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FOREWORD 


This  document  is  the  first  of  two  volumes  comprising  the:  final  report  on  a 
feasibility  study  of  techniques  for  the  control  of  re-entry  vehicles.  The  research 
was  sponsored  by  the  Flight  Dynamics  Laboratory,  RTD,  under  Task  No.  82181 
of  Project  No.  8225  and  Contract  No.  AF33(657)-7383  with  the  Minrieapolis- 
Honeywell  Regulator  Company,  Minneapolis,  Minnesota.  Lt.  E.  B.  Stear  was 
the  RTD  project  officer  for  this  program.  Work  at  Honeyv  ell  was  performed  by 
the  Systems  Techniques  group  of  the  Military  Products  Research  department  unde 
the  supervision  of  J.  T^  Van  Meter.  Project  personnel  were  :  G.  D.  Swanlund 
and  D,  K.  Scharmack,  principal  investigator s,  and  M.  D.  Ward,  R.  G.  Johnson, 
L.  D.  Dolid,  D.  L.  Lukes,  W.  C.  Marshall,  and  E..  R.  Rang. 

The  first  volume  of  this  report,  "Calculation  of  Optimal  Trajectories  and  . 
Synthesis,  of  Control  Functions",  presents  the  mathematical  models,  the 
development  of  theory  and  the  results  of  the  numerical  calculations  for  a  two- 
year  study  of  controlled  re-entry  vehicles.  The  digital  computer  programs 
developed  and  used  in  this  effort  are  documented  in  Volume  2,  "Computer  Pro¬ 
grams".  These  reports  are  designated  as  Honeywell  MPG  Documents  1535-TR2, 
Volumes  1  and  2. 

The  encouragement,  assistance,  and  direction  provided  by  Lt.  Edwin  B.  Stear 
is  gratefully  acknowledged. 

Some  of  the  digital  computation  was  conducted  at  the  RTD  computation  center. 


ABSTRACT 

Optimization  techniques  are  used  to  synthesize  the  control  programs  and  compute 
the  corresponding, flight  paths  for  controllable  re-entry  vehicles.  Linear  per¬ 
turbation  control  about  these  reference  trajectories  is  investigated. 

A  large  portion  of  the  theory  of  the  calculus  of  variations  is  modified  and  ex¬ 
tended  to  apply  to  this  problem.  Many  details  of  computational  techniques, 
necessary  in  the  adaptation  of  the  theory  to  large  scale  digital  calcuation,  are 
reported.  The  optimization  method,  which  evolves  is  a  modification  of  a  Newton- 
Raphson  iteration,  although  gradient  procedures  are  also  studied. 

The  -criteria  for  re-entry  trajectories  are  functionals  of  the  motion  related  to 
the  heating  of  the  v  ehicle.  It  is  found  that  these  criteria  are  relatively  insensi¬ 
tive  to  the  flight  path,  and  this  fact  leads  to  computational  problems  which  must 
be  handled  carefully.  The  paths  and.  control  authority  are  constrained  by 
reasonable  physical  requirements. 

The  linear  perturbational  control  is  found  by  requiring  the  integral  of  the  square 
of  the  control  deviation  to  be  minimum.  The  vehicle  position  is  the  object  of 
control  in  the  cases  studied,  it  is  found  that  this  form  of  control  can  be  used 
throughout  the  entire  re-entry  corridor.  Various  modifications  of  the  control 
gain  program  at  the  end  of  the  trajectory  are  considered. 

This  technical  documentary  report  has  been  reviewed  and  is  approved. 


Technical  Director 
Flight  Dynamics  Laboratory 
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SECTION  1 
INTRODUCTION 

COMPENDIUM 

Synthesis  of  control  for  a  .re-entry  vehicle  presents  many  difficulties  in  the 
form  of  norilinearities,  constraints  on  the  dynamical  state,  constraints  on 
the  control  and  questions  of  controllability  and  observability.  These  are  of 
such  importance  they  must  be  directly  incorporated  into  the  design  process 
and  cannot  be  considered  secondary,  as  is  done  in  most  engineering  methods. 
Further;  there  is  good  physical  justification  for  minimizing  certain  func¬ 
tionals  of  the  motion,  such  as  heat  input  and  acceleration  effects.  The  re¬ 
entry  problem,  therefore,  provides  an  excellent  proving  aground  for  the  new 
optimization  design  techniques. 

The  goal  of  this- study  is  to  demonstrate  the  efficacy  of  these  new  techniques 
by  actual  synthesis  and  simulation  of  the  mechanization  of  an  automatic  rer 
entry  control  scheme -developed;  through  their  use.  This  has  not  yet  been 
achieved.  However,  a  method  of  automatically  computing  optimal  re-entry 
trajectories  has  been  developed  and  demonstrated;  and  a  careful  study  of  a 
scheme  of  linear  control  about  a  reference  trajectory  has  been  conducted, 
serving  to, point  out  sensitivities:  and  difficulties  which  must  be  anticipated 
in  a  mechanization  attempt..  Many  details,  concerning  the  application  of 
the  theory  of  the  calculus  of  variations  to  this  problem  and  of  the  use  of  a 
large  digital  computer  in  trajectory  and  control  calculations,  have  been 
studied  and  are  reported  here. 


Manuscript  released  by  authors,  October  14,  1963,  for  publication  as  an 
FDL  Technical  Documentary  Report. 
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Model  £or  Re-entry  Studies 

The  model  chosen  for  these  re-entry  studies  is  two-dimensional  and  assumes 
a  spherical,  non-rotating  earth  and  an  exponential  atmosphere.  The  model 
differential  equations  are  a  great  simplification  of  the  actual  equations  of 
motion  but  do  represent  the  essential  phenomena  found  in  re-entry  flights. 

The  velocity  v,  flight  path  angle  y,  the  dimensionless  altitude  ?  -  h/R, 
and  the  great  circle  range  C  =  R  9,  are  chosen  as  the  components  of  the 
state  x.  The  control  is  represented  by  the  function  u  and  is  manifested  by 
authority  over  the  vehicle  aerodynamic  coefficients. 

In  this  notation,  the  .model  equations  of  motion  are 


dv  S  9  gn  siny 

If  ‘  -  tov  pv  CD’(u>  '  tftW 


d  y  S  v  cos  y  g  cos  y 

—  =  —  pv  C  (u)  +  -  -^—-2 

dt  2m  L  R(1  +?)  v(l  +  ?) 


(1.  1) 

d  ?  v 

—  =  —  sin  y 
dt  R 


dC  v 

—  =  -  cosy 

d  t  1  +  ? 


where  R  is  the  radius  of  the  earth,  S  and  m  the  vehicle  frontal  area  and 


mass,  respectively,  p  the  atmospheric  density,  and  Cj. 
aerodynamic  lift  and  drag  coefficients. 


and  Cj-j  are  the 


2 


For  convenience,  the  system  (1, 1)  is  usually  represented  by  the  vector 
equation 

x  =  f(x,  u)  .  (1.  2) 

Vehicles 


Three  vehicles  are  considered  in  these  studies.  The  first  is  a  high  lift-drag 
ratio  vehicle,  for  which  the  control  function  is  the  angle  of  attack  a  .  The 
aerodynamic:  coefficients,  based  on  flat  plate.:N  ewtorii  ah'  'flow ,  are 


CDO  +  CDL 


I  •  3 
I  sin 


sin  or  cos  a 


I  sin  at 


(1.3) 


in  which  C 


DO' 


CDL  an<3  ^LO  ai^e  constants.  The  other  two  vehicles  are 


low  lift-drag-ratio  Apollo-type  bodies.  One  has  a  variable  angle  of  attack,, 
for  which  the  coefficients  are 


CD  “  CDO  +  ^DL  cos  u  ’ 


CT  =  CT  ~  sin  u  ,  u  =  eat,  e  -  constant 

Lj  lu 


(1.  4) 


The  other  assumes  a  fixed  angle  of  attack,  but  uses  a  roll  angle  out  of 
the  plane  of  motion  as  the  control  function: 


CD  =  CDO 

CL  =  CLO  cos  ^  • 


(1.5) 


3 


The  criteria  for  those  optimal  trajectory  calculations  which  were  studied  are 
the  minimum  of  the  integral  of  heating  rate  deviation  from  a  prescribed  value 

T 

-jr-  J  (Q  -  q)2  dt  (U6| 

o 


and  the  minimum  of  total  heat,  represented  by 


6 


(T.  7) 


where  q  is  the  heating  rate  for  the  vehicle,  and  Q  is  a  given  constant; 
The  vehicle  described  by  Equations  (1.  3)  was  used  with  criterion  (1.  6); 
whereas  the  vehicle  represented  by  (1.4)  was  used  for  the  minimum  total 
heat  problem. 


Studies  Conducted 

A  considerable  body  of  theory,  based  on  the  calculus  of  variations,  was 
developed  to  facilitate  the  solution  Of  the  control  optimization  problem  on 
the  computer.  The  objective  of  these  studies  was  the  development  of  an 
automatic  optimization  method.  A  gradient  scheme  and  the  second-order 
Newton-Raphson  method,  using  the  first  and  second  variations,  are  des¬ 
cribed  in  detail  in  Section  II.  Also  presented  is>a  Newion-Raphson  scheme 
for  the  two -point  boundary  value  problem  associated  with  the  minimization 
problem.  Experiences  with  these  methods,  together  with  computer  results, 
are  given  in  Section  III. 


Linear  control  about  a  reference  trajectory,  found  experimentally  for  the 
roll-modulated  vehicle,  and  about  ah  extremal  trajectory  for  the  variable 
angle- of- attack  vehicle  is  derived  by  requiring  the  quadratic  integral  of 
control  deviation  Au, 

T 

f  (Au)?  dt  ,  (4.8) 

6 

to  be  minimum.  The  theory  arid  results  for  these  studies,,  along  with 
results  of  navigation  and  configuration  investigations,  are  given  in  Section  IV. 
A  detailed  error  analysis  and  a  sensitivity  study  of  the  controlled  systems 
also  are  outlined. 

This  volume  ends  with  conclusions  and  recommendations  in  Section  V.  The 
computer  programs  used  are  documented  nr  Volume  Two;.  ____ 


BACKGROUND  ~ 

Previous,  Studies  - 

Re-entry  into  earth  atmosphere  by  vehicles  launched  from  the  earth  or  from 
an  orbit  about  the  earth  is  a  problem  which  has  received  extensive  analysis 
during  the  past  several  years.  Allen,  Eggers,  Gazely,  Lees,  and  Chapman 
(References  1  through  7)  have  made  pioneering  contributions  to  the  under¬ 
standing  of  re-entry..  Their  important  work  has  been  augmented  by  a  multi¬ 
tude  of  studies,  most  significantly  those  of  Robinson  arid  Besonis  (Ref.  8), 
Wong  and  Sly e  (Ref,  9),  Young  and  Eggleston  (Ref.  10),  and  Luidens  (Ref.  11). 
A  comprehensive  list  of  re-entry  literature  would  be  extremely  extensive. 


The  primary  concern  of  the  early  studies  was  aerodynamic  heating  arid 
deceleration,  with  the  object  of  showing  the  feasibility  of  safe  descent 
through  the  atmosphere  from  orbital  initial  velocity.  The  earliest  studies 
concentrated  on  ballistic  vehicles.  Later,  the. effects  of  constant  lift- 
drag  ratios  were  evaluated.  More  recently,  careful  consideration  has 
been  given  to  modulated  aerodynamic  coefficients,  produced  by  varying 
the  angle  of  attack  according  to  some  stated  objective  (maintenance  of 
constant  deceleration  for  example);  References  12  through  17  are  examples 
of  this  work. 

Many  trajectories  have  been  calculated  for  a  wide  variety  of  initial 
velocities  and  flight-path  angles.  Attempts  to  generalize  solutions  through 
nonrdimensipnalization,  notably  the  work  of  Chapman  (Ref.  6  and  7),  and 
more  recently  Loh  (Ref.  18),.  have  been  particularly  rewarding. 


The  Safe  Flight  Corridor 


These  numerous  studies  resulted  in  clarification  of  the  concept  of  a  safe 
flight  corridor.  This,  corridor  is  usually  described  in  a  plot  of  altitude 
versus  velocity,  bounded  on  the  lower  side  by  heating  and/or  deceleration 
limits  and  on  the  upper  side  by  sustained  flight  at  maximum  lift  (negative 
lift  for  super -orbital  velocities  and  positive  lift  for  suborbital  velocities). 
Corridor  depth  is  used  to  evaluate  the  guidance  requirements  for  separa¬ 
tion  of  the  perigees  of  two  vacuum  conic  trajectories,  one  initiated  at 
the  upper  corridor  bound  and  the  other  at  the  lower;  this  is  described- by 
Chapman  (Ref.  7)  and  Luidens  (Ref.  11). 

As. a  result  of  this  work,  the  kinematics  and  aerodynamics  of  re-entry  are 
how  quite  well  understood. 
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Vehicle  Dynamics  and1  Control 

Recently,;  attention  has  been  given  to  vehicle  dynamics  and  control  during 
re-entry  (aside  from  a  few  early  studies  of  the  oscillatory  stability  of 
re-entry  bodies).  Several  relatively  simple  range  and  cross-range 
schemes  have  been  proposed  for  landing-point  control.  References  16  and 
17  are  examples  of  this  kind  of  investigation,  based  on  currently  established 
concept s-andt equipment  mechanization. 

Gredit  for  the  use  of  optimization  techniques  in  determining  re-entry 
trajectories  seems  due  Bryson,  et  al  (Ref.  19).  This  paper  pointed  out  the 
possibility  of  using  more  sophisticated  control  over  the  vehicle  (to  minimize 
functionals  of  the  motion)  while  still  satisfying  given  end-point  conditions. 
Linear  control  techniques  were,  again,  first  published  by  Bryson  (Refs.  20 
and  21).  These  methods  provide  the  vernier  adjustments  to  the  control 
function  made  necessary  by  model  variations,  such  as  atmospheric  density 
deviations  and  winds. 

The  control  scheme  envisioned  in  this  report  .makes  use  of  both  optimal 
trajectories  and  linear  control.  Since  Bryson's  and  similar  calculations 
are  not  automatic  but  require  the  intervention  of  the  person  making  the 
calculation  to  produce  convergence  of  the  iterations,  these  methods  are 
not  suitable,  without  modification,  for  on-board  mechanization.  The 
studies  reported  here  were  aimed  at  finding  a  calculation  which  would  be 
automatic,  converge  rapidly  and  ensure  an  actual  minimum  for  the 
functional. 
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}  SECTION  II 

‘t  OPTIMAL  TRAJECTORY  CALCULATIONS  -  THEORY 

% 


1  PROBLEM  STATEMENT  (A) 

% 

<■  _  A 

*  The  problem  considered  here  is  a  special  form  of  the  problem  of  Bolza  as 

<■  formulated  by  Bliss  (Ref.  22)  and  extended  to  include  inequality  constraints 

by  Valentine  (Ref.  23).  It  is:  Find  that  path  which  minimizes  the  function 

T 

J  -  (T,  x  (T))  +  J  iQ  (x,  u)  dr  (2.1) 

o 

subject  to  differential  equations  of  the  form 

x  -  f  (x,  u),  x  (  0  )  =  xq  ,  (2.,2) 

inequality  constraints 

G  (x,  u)  a  o  ,  (2,  3) 

and  terminal  surface  equations 

^(T,  x  (T))  =  0  .  (2,4) 

In  the  above-,  x  and  u  are  n  and  m  dimensional  state  and  control  column 
vectors,,  and  (*)  represents  differentiation  with  respect  to  the  independent 
variable  t.  The  dimension  of  vector  Equation  (2.  3)  is  q,  and  that,  of 
(2. 4)  is  r,  where  r  £  n  +  1. 
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(2.5) 


Following  Valentine  (Ref.  23),  Equations  (2.  3)  are  rewritten  in  the  form 

a2  =  G  (x,  u) 

•  2  *2 

where  the  components  of  the  vector  a  are  o'.  ,  j  =  1,  q.  The 

'I* 

•o*s  are  slack  variables, introduced  to  complete  the  set  of  differential 

v  V  •  .  < 

equations  for  the  problem  of  Bolza.  No  more  than  m  of  the  a*s  may  be 
zero  at  .any  point  on  the  path. 


It  is  assumed  that  f  (x,  u)  and  the  differential  Equations  (2.  2)  and  (2.  5) 
have  continuous  partial  derivatives  of  at  least  third  order  in  all  variables, 
in  an  open  region  S.  about  the  minimizing  path.  Furthermore,  the  matrix 
made-up  of  the  partial  derivatives  of  the  differential  equations,  with  respect 
to  all' derivatives  and  the  control  functions,  must  have  rank  .n  +  q  at  each 
point  of  the  minimizing  path.  This  ensures  that  the  differential  equations 
are  independent.  The  matrix  has  the  form 


”  I  V  Z 


’  Z’  VuG  -2 


(2.6) 


where  I  is  the  nxn  identity  matrix,  Z  is  an  nxq  zero  matrix,  Z’  is- 

the  transpose  of  Z,  and  2.  is  a  qxq;  diagonal  matrix  whose  elements 

are  2a... 

J 


The  solution  x  (t)  of  the  differential  Equations  (2.  2)  and  (2.  5)  is  supposed 
to  be  continuous,  with,  at  least,,  absolutely  continuous  first  derivatives. 

In  some  instances,  it  will  be  possible  to  consider  corners,  i.  e. ,  points 
at  which  the  derivatives  are  discontinuous.  The  control  functions  are 
treated  as  derivatives  in  (2.  6)  so  they  can  be  discontinuous,  as  in  the 
bang-bang  problem.  Potential  corners  are  those  points  at  which  inequality 
constraints  change  from  greater  than  to  equality  states  or  vice  versa.. 
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They,  may  also  be  defined,  by  switching  points,  as  in.the  bang-bang  problem 
where  the  equality  sign  holds  in  the  constraint  relation.  A  subarc  is  de¬ 
fined  as  that  part  of  the  path  between  corners  or  potential  corners.. 


Finally,;  the  functions  g  (T,  x(T)).  and  (2.  4)  are  assumed  to  have  continuous 
partial  derivatives  of  at  least  third  order  in  an  open  set  S of  points 
(T,  x  (T)),  and  the  matrix 


-3.^ 

3  T  b  Xj. 


3^ 

3xn 


Y.l// 


(2.  7) 


0% 

oT  ’  3x1  * 


at- 


3x 


n 


is  assumed  to  have  rank  r.  An  admissible  arc  is  defined  to  be.  a  path  having 
all  of' its  elements  (t,  x,  x,  u)  and  (T,  x  (T) )  lying  in  S.j  and  S2> 
respectively. 


Differential  equations  and  inequality  constraints  containing  t,  explicitly, 
can  easily  be  brought  to  the  form  of  Equations  (2,  2)  and  (2.  5)s  The 
independent  variable  is  changed  from  t  to  s  by  adding  the  differential  equation 


dt 

ds 


and  noting  that 


dx 

dx 

ds 

dx 

dt 

ds 

dt 

ds 

(2.8) 


(2.9) 
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Then  t  becomes  a  dependent  variable,  and  the  resulting  n  +  q  +  1  differential 
equations  are  in  the  desired  form.  Observe  that  the  order  of  the  system;  har, 
been  increased  by  one. 

NECESSARY  CONDITIONS  (B) 

The  necessary  conditions  for  this  problem  are  stated  here  without  proof.  The 
reader  is  referred  to  References  22,  23,  and  24,  for  their  derivation.  The 
first  necessary  condition  is 

The  Multiplier  Rule 

An  admissible  arc  E,  defined  on  an  interval  [  0  ,  T]  is  said  to  satisfy  the 
multiplier  rule  if  there  exist  constants  pQ  =  1,  e’  =  [  e^,  .  ,  e^j-  not 

all  zero,  and  a  function 

F(t,  x,  u,  x,  p,  n,  b)  =  fQ  +  p*  (f  -  x)  +  m- *  (G  -  a2)  (2.  10) 

with  multipliers  p'(t)  =  ( p.  (t),  ...  ,  pn<t)],  continuous  on  [  0  ,  T  ]  and 

lx*  (t)  =  [(J..  (t),  ...  ,  m-  (t)  ]  continuous  on  [0  ,  T]  except  possibly  at 

corners  of  E  where  unique  right  and  left  hand  limits  exist,  and  satisfying 
the  Euler-Lagrange  equations 

-p*  =  y  f  +  p'  M  +  [i'V  G  (2.  11) 

^  X  O  x  r  X 

0  =  V  f  +  p»  V  f  +  v  G  (2.  12) 

0  =  (i.  G.,  j  *  1,  ...  ;  q  .(2.13)- 

J  J 

where 

0  *  n  ,  <2. 14} 
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v^?> 


arid  differential  equations 


x  =  f; 


a  =  G  s?  0 


x  (0)  =  x 


(2.  15) 
(2.16) 


along  E,  and  furthermore,  such  that  the  equations 


f  +  p'f  4- 
=o 


g 


dT 


+  e* 


b\p, 
S  T 


•dT 


(2.  17) 


+  (Vxg  +  e!V^  -  p‘) 


dx  (T)  =0,  ip  =  0 


hold' for  every  choice  of  the  differentials  dT,  dx(T).  The  multipliers 
PQ,  P,  l1  do  riot  vanish  simultaneously  at  any  point  of  the- interval  [  0,,  r 
for  an  arc  E  satisfying  the  multiplier  rule.  Furthermore,  the  function 


H  s  f  +  p»f 
o 


is  a  constant  on  E.  Every  minimizing  arc  E  for  the  given  forth  of  the 
problem  of  Bolza  must  satisfy  the  multiplier  rule.  In  the  above,  the 
is  the  first  column  of  the  matrix  (2;  7). 

T 


vector 


dip 

bT 


(2.  18) 


The  vanishing  of  the  coefficients  of  dT  and  dx  (T)  in  equation  (2.  17) 
constitutes  the  transversality  condition  on  the  arc  E,  Thus, 


H  =  <fQ  +  p«f)t  =  (fQ  +  ,p'f). 


g 


3  T 


+  e’ 


b  ip 


bT 


(2.  19) 


P‘(T)  =  (V  g  +  ef  V  ^)-i  . 

^  !  T 


(2.20) 
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The  Hamiltonian,  Equation  (2. 19),  is  often  called  a  first  integral  for  this 
problem  It  has  further  significance  in  the  second  necessary  condition: 


l 


The  Necessary  Condition  of  Weierstrass  (The  Minimum  Principle) 

An  admissible  arc  E  satisfying  the  multiplier  rule  with  multipliers 
PQ  1»  P*  p#  is  said  to  satisfy  the  necessary  condition  of 
Weierstrass  with  these  multipliers  if  the  condition 


H  (t,  x,  p,  u)  £  H  (t,  x,  p,  U) 

is  valid  at  every  element  (t,  x-,,  x,;  u)  of  E  for  ail  admissible  points; 
«• 

(t,  x,  X,  IJ)  ^  it,  x,,  x,  u)  satisfying  the  Equations  (2.  2)  and  (2.3). 
Every  minimizing  arc  E  for  the  given  form  of  the  problem  of  Bolza 
must  satisfy  this  condition. 

A  consequence  of  the  Weierstrass  condition  is: 


The  Necessary  Condition  of  Clebsch 

r  "  "  " 

At  each  point  of  E  let  G  be  a  vector  whose  components  are  those 
r  components  of  G  which  vanish,  and  {x  the  corresponding  multipliers. 

Then  the  inequality 

L  n*  (H  +  pG)  tt  *0 

r 

L 


I 
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(2.  21) 


(2.22) 


must  be  satisfied1  for  every  vector  fr  ^  0  where  tt*  -  ,  *  *  ■>  •*-  nm-i 

and  satisfies  the  equations 

V  £  tt  =  0.  (2,23) 

u 

Every  minimizing  arc  for  the  given  form  of  the  problem  of  Bolza  must 
satisfy  this’  Condition. 


EXTREMALS  AND  THE  EQUIVALENT  MINIMIZATION  PROBLEM  (G) 

The  paths  considered  here  are  supposed  to  satisfy  the  Equations  (2.  11) 
through  (2.  16),  as  well  as  the  Weierstrass  and  Clebsch  necessary 
conditions*  Such  paths  will  be  called  extremals.  The  paths  are  further 
supposed  to  consist  of  a  finite  number  of  subarcs*  On  each  subarc  a 
given  subset  of  the  inequality  constraints  are  equality  constraints,  and 
all  the  rest  are  greater  than  zero,  except  possibly  at  a  finite  number  of 
points.  It  is  shown  in  Appendix  A  that  such  paths  can  be  generated  by 
integrating  sets,  of  differential  equations  of  the  form 


Sc  =  f  (x,  p)  „ 

x (0)  =  x 
-0 

(2.  24) 

-P  =  v^F(x,  p) 

(2.25) 

These  differential  equations  change  from  subarc  to  subarc,  but  the 
solutions  are  continuous  at  the  junction  of  two  subarcs.  It  is  further 
shown  in  Appendix  A  that  the  solutions  can  be  represented  as  functions 
of  the  initial  conditions  On  the  multipliers  p(0)  =  pQ,  and  the  terminal 
value  of  the  independent  variable,  T;  i.  e. ,  a,  solution  is  defined  by 
specifying  a  particular  set  (T,  pQ)..  Continuous  partial  derivatives  of  at 


*  Junction  points  are  potential  corner  points  and  are  also  referred  to 
as  breakpoints. 
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least  second  order  in  these  variables  exist.  The  optimization  problem  can 
then  be  reformulated  as  :  Minimize  the  function 

J  =  J  (y)  (2.26) 

in  the  variables  y  =  (T,  pQ),  subject  to  the  constraint  equations 

4f(y)  =  0.  (2.27) 

This  is  the  equivalent  minimization  problem.  Necessary  and  sufficient 
conditions  for  a  relative  minimum  for  this  problem  are  well  knowjn. 

First,  form  the  function 


F  (y,  X)  =  J  (y)  +  X*  s//(y)  ,  (2.  28) 

where  X  is  an  r  dimensional  vector  of  Lagrange  undetermined  multipliers. 
Then  assuming  that  y  is  defined  Over  an>  open  region  R,  that  the  necessary 
partial  derivatives  exist  in  the  neighborhood  of  the  critical  point  y  =  c, 
and  that  the  critical  point  is  normal  (Ref.  22,  pp,  210-213  .),  the  necessary 
and  sufficient  conditions  may  be  written: 

If  a  critical  point  y  =  c  has  a  set  of  multipliers  X  for  which  the 
function  F  satisfies  the  conditions 

V  F  (c)  =  0  (2.29) 

y 


and 

Ay»  V2  F  (c)  Ay  >  0  (2.30) 

y 

for  all  vectors  Ay  satisfying 

V  Ay  =  0  ,  Ay  ^  0  ,  <2. 31) 

J 

then  c  is  a  minimizing  point. 
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It  should  be  remarked  that  the  Equations  (2.  29)  and  (2,  27)  are  n  +  r  +  1 
equations  in  the  n  +  r  +  1  variables  y  and  X.  They  may  be  solved,  in 
principle  at  least,  for  the  critical  value  y  =  c  and  the  corresponding 
set  of  multipliers,  if  the  matrix 

V2  F  (V  tf/)» 

y  yr  ; 

(2.32) 

V  \f/  Z 

y 

is  non- singular. 

The  sufficiency  condition.  Equations  (2.  30)  and  (2.  31)  is  rather  awkward  to 
test  numerically.  However,  a  more  convenient  test,,  showing  that  a  matrix 
is  positive-definite,  is  easily  derived.  There  are  two  ways  in  which  this 
can  be  done. 

The  most  obvious  method  is  to  solve  (2,  31)  for  a  set  of  r  dependent  components 
of  Ay  in  terms  of  the  (n  +  1  -  r)  others,  and  to  substitute  back  into  (2.  30). 

The  resulting  (n  +  ,1  -  r)  x  (n  +  1  -  r)  matrix  must  be  positive-definite.  The 
other  method  is  to  complete  the  set  of  constraint  Equations  (2.  27)  with 
(n  +  1  -r)  equations  of  the  form 

(y)  =  ^  (c)  +  z  (2.33) 

chosen  such  that  the  matrix 


is  non- singular. 


(2.  34) 
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Then  the  equations 


j  ? 


;  -o- 

i  az 


r  <2. 35) 


may  be  thought  of  as  a  transfer  ma ti bn  r elafcin  g  ar  bitrary "  A„z  to  values 
of  Ay  satisfying  the  constraint  Equations  (2,  31).  Substitution  of  5); 
into  the  sufficiency  condition  (2.  30),  followed  by  the  indicated 
multiplication  then  gives  the  result  ..  • 

Az1  Aa  Az  >  6-.,  (2.36) 

4  " 


Since  Az  is  arbitrary,  the  (n  +  1  r  r)  x  (n  +  1  -  r)  matrix  A^ 
must  be  positive-definite. 


COMPUTATIONAL  METHODS  FOR  FINDING  A  RELATIVE  MINIMUM-(D) 


Three  of  the  many  numerical  schemes  for  finding  a  minimum  are  discussed 
here.  These  are  the  Newton -Raphson  method,  a  modification  of  the  Newton- 
lethod  and  the  method  of  steepest  descent.  Reference  31  gives  a 
other  methods.  It  i's  assumed  that  the  required  partial  derivatives 
are  available.  Methods  for  obtaining  them  are  given  in  subsections  E  and  F. 

Consider  first  the  case  in  which  r  ,  the  dimension  of  terminal  surface 
Equations  (2.  4)  is  less  than  (n  +  1).  Expansion  of  Equations  (2.  26)  and  (2. ,27) 
through  second-order  terms  about  a  point  yQ  gives 

J(y)  =  J.(y0)  +•  VyJ<y0)  Ay  +  y  Ay»  j(yQ)  Ay  +  ...  (2.37) 

Wy)  =  ^(y0)  +  vy^(y0)Ay  +  YAy'  Vy  ^(yo)  Ay  +*"  *  (2*38) 

where 

^y  =  y  -  y0* 
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fe  thj  Ne^tori-Haphson  scheme;  the  point  yQ  is:  assumed  to  be  so  close  to  the 
minimising;  point  that  (2.  37)  and  (2.  38)  adequately  approximate  behavior  at 
the  minim iziiig ■point,  The  function  F  of  Equation  (2.  28)  is  formed  with  these 
expressions  and  differentiated  with  respect  to  y,.  the  assumed  minirnir ing 
point.  Setting  the  result  to  zero  gives 

7y  P(y,  X.)  =  0  =  '?  P(yo>  X)  -f  V2  F(yo,  X)  Ay,  (2.39) 

The  iterative  solution  of  Equations  (2.  39)  and  (2.  38)  (with  ift( y)  =  0)  for 
Ay  and  \  then  constitutes  the  Newton-Raphson  method,  provided  the  matrix 

F(y0.  X)  <2-40> 

satisfies  the  sufficiency  condition.  Otherwise,,  the  solution  would  be  driven 
toward  a  saddle  point  or,  worse,  a  maximizing  point. 

If  the  sufficiency  condition  is  not  met,  some  other  iterative  scheme,,  such  as 
steepest  descent,  should  be  used.  This  situation  will  generally  arise  if  the 
point  yQ  is  too  far  from  the  minimizing  point  for  the  second  order  expansions 
(2.  37)  and  (2.  38)  to  be  accurate.  In  steepest  descent,  the  point  y  is  assumed 
to  be  closer  to  the  minimizing  value  than  the  present  estimate  y  ,  The  ex¬ 
pansions  (2.  37)  and  (2.  38)  are  truncated  after  first-order  terms  and  written 
in  the  form 

r  dj 

di // 


:  = 

J(y) 

-  •%/ 

vy  J(*o> . 

i 

-  ^(y0) 

Ay  =  9  ( yQ)  Ay  (2.41) 
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A  change  Ay  in  the  gradient,  or  steepest  descent,  direction  is  made  by  the 
choice 

Ay  =  9 1  (yQ)  K?  (2.  42) 

where  K  is  an  (r  +  1)  x  1  constant  vector  to  be  determined.  Substitution 
of  (2. 42)  into  (2,  41)  to  determine  K,  and  insertion  of  the  result  back  into 
(2.  42)  gives 


;  dJ 

;  dp 


(2.  43) 


that  is,  the  desired  change  in  yQ  is  .given,  in  terms  of  the  G  (yQ)  matrix, 
arid  specified  small  changes  dJ  arid  dp  .  The  latter  are  chosen  to  make 
J  and  p  smaller  at  the  completion  of  the  next  step.  Reference  (19) 
provides  a;  more  complete  description  of  the  use  of  the  method.  It  should 
be  pointed  out  that  the  method  blows  up  as  a  critical  point  is  approached, 
since  the  (r  +  1)  x  (n  +  1)  matrix  6  is  necessarily  of  rank  r  at  such  a 
point.  Otherwise  the  gradient  vector  could  not  be  the  zero  vector.  The 
inverse  matrix  of  Equation  (2.  43)  does  not  exist  under  this  condition. 

Hopefully,  the  critical  point  is  a. 'minimizing- -point,  so  the  computational 
equations  can  be  switched  to  those  of  the  New.ton-Raphsdn  method  as  the 
point  is  approached. 

Now,,  suppose  r  =  (n  +  1).  In  this  case,  the  (n  +  1)  Equations  (2.  38) 
completely  specify  the  (fixed)  end-point  and,  hence.  Ay,  without  reference 
to  the  minimization  criterion  (2,  37);  Truncation  of  (2.38)  after  first-order 
terms  and  rearrangement  then  gives 

Ay  =  [  Vyi^yJ]"1  Wy>  -  i//(yQ))  .  (2.44) 
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provided  the  problem  is  normal.  If  the  choice 

tf%)  -  ^(yb)  =  “  G  iJ/(yo) 

4 

is  made,,  a  gradient-like  iterative  equation  results:. 

Ay  «  e  I  Vy  ti/(yo) )  -1  ^(yo)  ,  0  <  C  £  1 

Use  of  an  equation  of  this  form,  to  find  the  minimizing  path  will  be 
termed  the  "modified  Newton- Raphson  method.  " 

The  scalar  constant  C  is  the  only  experience  factor  required  for  this 
scheme.  If  C  =  1,  Equation  (2.  46)  is  recognized  as  the  ordinary 

Newton-Raphson  method  for  finding  a  root.  Usually,  the  initial  guess 
of  yQ  is  so  far  from  the  minimizing  .point  that  some  smaller  value  for 
C  will  be  required;  hence,  the  bounds  0  <  C  4  1.  As  the  solution 
approaches  the  critical  value,  it  is  normally  found  that  acceptable  values 
of  C  become  larger  and  larger.  It  is  possible  to*  automate  the  proper 
selection  of  C  on  the  computer. 

There  is  no  guarantee  that  the  method  will  converge  to  a  minimizing 
solution.  However,  the  second  variation  test  described  in  Appendix  B 
can  be  used  to  determine  whether  or  not  a  solution  does  minimize. 

The  modified  Newton-Raphson  method  can  also  be  used  for  the  class 
of  control  optimization  problems  for  which  the  terminal  surface 
Equations  (2.4)  assume  either  the  form 

i//.  =  x.  (T)  -  X.  =  0  ,  i  =  1,  ....  r-  <  (n+  1.) 


(2,  45) 


(2.  46) 


(2.47) 
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or  the  form 

■>  ' 

'k  =  X  (T)  r  X  =  0  ,  j  -  1,  ...  ,  (r  -  i)  <  n  (2.48) 

J'  r}  J 

\pr=  T  -  D  =0,  (2.49) 

ije 

after  a  suitable  renumbering  of  the  state  coordinates;  In  these 

equations,  X.  ,  X.  and  D  are  specified  constants.  A  large 
*  1 

number  of  optimization  problems  fall  into  this  category.  The 
-advantage  of  using  the^  modified;  New.tOn-Raphson>  method. over  the 
method  of  steepest  descent  is  faster  convergence. 

When  the  terminal  conditions  (2.  47)  are  substituted  into  the  transversality 
conditions  (2. 19)  and  (2.  2.0)  „  there  results 


H  <V  P6J  ' 

(f 

O 

4  P'O  i  =  0  ■  0 

<2.  30) 

-p.  (T) 

er 

i  =  1,  .  ,  r 

(2.  51) 

p.(T) 

0,. 

j  =  r  +  1,  ...  n 

(2;  52) 

Since  the  Hamiltonian  is  constant  along  the  entire  path  it  is  evaluated 
in  (2.  50)  at  the  instant  t  =  0,  which  allows  simpler  computation. 
Equation  (2.  51)  specify  the  constants  e^  and  hence  give  no  new  infor¬ 
mation.  When  the  solutions  are  considered  to  be  functions  of  the 
vector  y,  Equations  (2.  47),  (2.  50)  and  (2.  52)  become  the  set 


Xf  (y)  -  X.  =  0  , 

i  =  1,  .  . .  ,  r 

H  (y)  =  0, 

(2.53) 

Pj  (y)  =  0  , 

j  =  r  +  1,  ...  ,  n 

*The  function  g  of  Equation  (2.  1)  is  taken  to  be  zero  here  for  convenience 
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This  is  the  two-point  boundary  value  problem,,  based  on  the  necessary 
conditions  of  the  minimization  problem..  When  the  set  (2.  53)  is  expanded 
about  the  point,  yQ  and  rearranged  as  in  th'e  development  of  Equation  (2.  46) 
there  results 


Ay  -  -OC!'V0) 


*<y0)  -  x 


K(ya) 


0  <  C  s  1„ 


(2.  54) 


where  (x  -  X)  and  p  are  the  vector  forms  of  x,.  and  p..  in  (2.  53). 
Let  f (0)  be  the  right-hand  side  of  Equation  (2.  2)  evaluated  at  t  =  0. 
Then  the  matrix  n  is 


n(y0> 


v  x(y  > 
y  'o 

0  ,  f*  (0)^ 


(2.  55) 


Again,  there  is  no  guarantee  that  the  iterative  Equation  (2.  54)  will 
converge  to  a  minimizing  solution.  However,  the  sufficiency  condition 
of  Subsection  C  may  be  used  to  determine  if  an  optimal  trajectory  has 
been  obtained. 


When/the  terminal  value  of  the  independent  variable  is  specified  as  in 
Equations  (2.  48)  and  (2.  49),  a  reduction  in  the  size  of  the  system  is 
found.  The  boundary  value  problem  for  the  terminal  surface  becomes 
that  of  solving  the  system 

x.(p  )  -  X.  =  0  ,  j  =  1,  . . .  ,  (r  -  1) 

J  ^  J 

Pk<PQ)  =  0  ,  k  =  r,  . . .  ,  n 


(2.  56) 
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for  the  minimizing  pQ  .  The  iterative  equation  is 


r  * 

x  x  <p0)  -  x 


Apo  r  -  G  w"  <po) 


0  <  C  *  1  ,  (2.  57) 


p  (p0) 


where  (x  ■?  X)  and  p  are  the  vector  forms  of  the  left-hand  side  of 
system  (2.  56),  and 


f  V  &(p  )  - 

Po  0 


w(po)'  = 


(2.58) 


V ..  p  (p  ) 

p0.  0 


THE  PARTIAL  DERIVATIVES  I  (E) 


When  Equations  (2.  26)  and  (2.  27)  are  differentiated  with  respect  to  the 
variable  y  =  (T,  pQ),  there  results 


dg(T,  x  (T), ) 

j  =  - - - - -  +  V  g  (T,  x  (T))  x(T)  +  f  (T) 

5  T  x  ° 


(2.59) 


V  J  =  V.  g(T,  x(T))  V  x  (T)  +  V  x  .  (T)  (2. 60) 

X  6  *  n  D_n+1 


o  \p(T,  x  (T)) 

ip  =  -  +  7  i//(T,  x  (T))  x  (T) 

b  T  X 


(2.  61) 


V  \b  =  V  MT,  x  T  V  x  (T) 
p  r  x T  p 

“o 


(2.62) 
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where 


wn+  1 


<u 


f  dT 

O  - 


y 


(2.  63) 


the  integral  term  of  Equation  (2.  1), has  been  introduced  as  a  new  variable 
(see  Appendix  A). 

♦ 

The  symbol  J  in  Equation  (2;  59)  represents  the  partial  derivative  Of 
Equation  (2;  26)  with-  respect  to  the  variable  T.  This  variable  appears 
both  explicitly  and1  implicitly  (in  the  variables  x),  as  indicated  on.  the  right- 
hand  side  of  (2.  59).  A  similar  interpretation  holds  for  Equation  (2.  61). 

Since  all  terms  of  Equations  (2.  59)  and  (2,  61)  are  available  after  a  path 

#  « 

has  been  found,  J  and  i//  are  readily  evaluated.  The  other  partials 

may  be  estimated  or  they  may  be  computed  if  7  x  (T)  and  7  x  (T) 

P-.  ‘Pa  n  +  1 

are  known.  u 


Consider,  first,  the  estimation  of  the  partial  derivatives.  Since  the 

solutions  are  continuous  functions  of  p  ,  the  function  J  may  be  expanded 

in  a  Taylor’s  series  expansion  through  second-order  terms  in  one  of  the 

variables  p,-  ; 

o 


Jlp.  4  5p.  ) 
o  o 


J  + 
o 


d  J 


d 


2 


J 

~2 

o 


<6p,  )2 

o 


(2.  64> 


The  left-hand  side  notation  of  (2.  64)  indicates  that  only  one  component  of 

pQ  has  been  changed  by  the  small  positive  amount  6p. 

1o 
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Note  that  the  terminal  value  of  the  independent  variable  T  is  also  held 

m  i2j 

constant.  The  quantities  J  ,  —  and;  v-,  -  y  correspond  to  the 

•  iq  ' 

extremal  generated  with  the  initial  conditions  pQ  .  Substitution  of 

-6p.  into  (?,;'64);  results  in 
'o 


"  5p.  )  =  JQ  - 


(5P.  ) 

o 


(2.65) 


Subtraction  of  (2,  65)  from  (2.  64)  and  rearrangement  gives,  the  desired 
result 


j(p  +  6p.  )  -  J(p.  -  5p;  ) 
o  _  o _ o  .  6 

2  6p. 

lb 


(2.  66) 


that  is,  the  partial  derivative  is  obtained  from  two  solutions  with  slightly 

different  initial  conditions*  in  one  of  the  components  of  pQ  .  A  complete 

set  of  first  partials  is  obtained  from  2n  perturbed  solutions.  Since 

values  of  can  be  determined  for  these  per-  jr.bed  solutions,  the  obvious 

changes  in  Equation  (2.  66)  may  be  made  to  calculate  V  \p  . 

po 

If  the  Newton-Raphson  method  is  used,  second  partials  must  also^be 
estimated.  The  form  of  the  matrix  of  second  partial  derivatives  for  the 
function  J  is 
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Of,  these,  only  J  can  be  calculated  explicitly.  The  computational  equation, 
from  differentiation  of  Equation  (2, '59),4s 


3  g(T,  x  (T))  3 

j  =  - - +  2  - (7vg(T,  x  (T))  x  (T)  +  7  g(T,  x  (T))  x  (T) 

3  f Z  3  T  x 


<2.  67) 


+  x*  (T)  vJgiT,.  x  (T))  x  (T)  +  fQ  (T)  , 


where 


f  =  v  f  x  +  7  f  (V  u  x  +  7  up), 
o  x  o  u  o  x  p  r 


(2.  68) 


and  where  7  t  and  7  f  are  partials  with  respect  to  explicit  appearance 
X  o  U  O  -  * 

of  x  and  u  in,  f  Similar  equations  may  be  written  for  x  =  f  .  The 

o  , 

components  of  7^  J  are  computed  from  application  of  (2.  66)  to  the 
function'  J  ... 


The  diagonal  elements  of  7  J  ar.e  estimated  from 

po 


32J  J(p.  +  6p.  )  +  J(p.  -  8p.  )  -  23 


(2.69) 


■<6Pi  Y 

o 


Equation  (2.  69)  is  derived  by  adding  Equation  (2.  64).  to  (2,  65)  and 
rearranging  the  result. 


The  off-diagonal  terms  of  V"  j  are.  estimated; from 

Po 


j<p>  +  '6p,  p.  +  dp.  )  +  J  -  J(p;  +  6p.  )  -  J(’p.  +  6p.  ) 

**  V 


•Bp.  B.p^ 


O  '  O; 


o  Jo 


$P-  6p. 

Jo 


(2.  70} 


where  J(p.  +  6p.  ,  p.  +  6p.  )  is  the  solution  which  results  when 

o  1o  *^0  ^o 

the  elements  p.  and  p.  are  perturbed  by  the  small  positive  amounts 
o  •'o 

6p.  and  6p.  ,,  respectively. 

1q  ^o 

Equation  (2.  70)  was  derived  by  expanding  J  in  a  Taylor’S' series  expansion 
through  second-order  terms  in  the  variables  p.  ,  p,  about  the  solution 

Ao  JO 

for  pQ,.  Substitution  of  (2.  66)  and  (2.  69)  into  the  result,  followed  by  some 
rearrangement,  then  gives  (2.  70). 

Equations  similar  to  (2.  70)  can  be  found  for  either  or  both  of  6p.  and>  6p. 
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considered  negative.  It  is  recommended  that  cross  partials  be  computed 

for  the  perturbations  6p.  ,  5p.,  taken  in  the  direction  of  expected  change 

“O  ^o 

in  the  solution.  This  ~ an  be  found  from  the  direction  of  steepest  descent.. 
The  partials  (2.  70)  require  an  additional  —  % — —  solutions.  The  grand 

total  of  solutions  for  a  Newton-Raphson  step  rises  to  — - -  ,  not  counting 

the  unperturbed  solution 

Finally,  second  derivatives  of  i/a  are  computed  by  substituting  components 
of  i//  for  J  in  the  above  equations,  with  the  exception  that 


* l  =  — r  +  2 


Vl]  x  +  vx, //.  x  +  x1  v2  gi 
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There  are  several  disadvantages  associated  with  estimated  partial  derivatives. 
The  biggest  problem  involves  the  integration  accuracy.  Since  the  partial 
derivatives  are  approximated  by  differencing  integrated  solutions,  the  random 
error  generated  by  the  integration  algorithm  may  result  in  bad  predictions 
of  the  minimizing  point.  Closely  related  to  this  problem  is  the  one;  of  picking 
the  perturbations,  6  pQ.  If  these  are  too  small,  the  partial  derivatives  will 
reflect  the  integration  error  alone.,.  If  they  are  too  large,,  the  partial  deriva¬ 
tives  will  not  reflect  the  nature  of.  the  surface  in  the  vicinity  of  the  point 
(T,  p  ).  Intermediate  values  must  be  chosen  by  experience:  Even  then,. 

O  '  *  *"  -  v  ‘  V'  -  - 

there  are  errors  in  the  estimating  equations,  since  these  corresporid  td^dii- - 
ferehtiated  (quadratic)  three-point  fits.  It  generally  takes  longer  to  estimate 
the  partials  on  the  computer  than  to  compute  them  explicitly  (see  Subsection  E). 
The  total  number  of  equations  integrated  per  step  is  n(n  +  1)  (n  +  2)  for  esti¬ 
mated  partials,  whereas  it  is  2n  (p  +  1)  for  explicit  computation..  It  took  about 
three  times  longer  to  estimate  the- partials  than  it  did  to  compute  them  explicitly 
for  a  particular  problem  ,(n.  =  3 ). 

In  spite  of  these  difficulties  it  may  be  wise  to,  estimate  the  partial  derivatives 
for  some  problems  rather  than  to  calculate  them  explicitly.  The  set  of 
differential  equations  for  estimated  partials  is  much  smaller  and  easier  to  set 
up  than  the  corresponding  set  for  the  explicit  computation  of  partials.  Even 
if  the  partials  are  to  be  computed  explicitly  it  is  helpful,  for  program 
checkout  purposes,  to  have  estimates  of  the  partial  derivatives. 

The  value  of  the  Hamiltonian  is  theoretically  constant,  so  it  may  be  used  to 
detect  large  truncation  and  round-off  errors. 


THE  PARTIAL  DERIVATIVES  II  (F) 

Now  consider  the  explicit  computation  of  partial  derivatives.  For  this 
purpose  it  is  convenient  to  shift  to  the  Hamiltonian  formulation  of  the 
problem.  It  may  readily  be  verified  that  the  equations 
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r 


m”f  Mi 


X  a  (x,  U,  p,  fi), 

p  =  (x,  u,  p,  JJL ) 

0  =  j  (x,  u,  p,  ja  ) 

Hj  =  H(x,  u,  p)  +  |iG(x,  u)  ,, 


(2.72) 


(2.  73) 


(2.  74) 


(2.  75) 


where  H  is  defined  by  Equation  (2.  18),  are  identical  with  Equations  (2.  15), 
(2,  11)  and  (2.  12),  respectively.  According  to  the  arguments  of  Appendix  A, 
each  subarc  of  a  path  has  its  own  set  of  reduced  differential  equations  of  the 
extremals  of  the  form 


x  =  V!  H1  (x,  p) 


(2.  76) 


-p  =  ^Ha  (x,  p). 


(2.  77). 


where,  of  course. 


u  =  u  (x,  p) 


(2.  78) 


li  =  p(x,  p)  , 


(2.  79) 


Now,  in  the  solutions  of  Equations  (2.  76)  and  (2.  77),  x,  p,  are  functions 
of  the  independent  variable  t  and  the  initial  conditions  x^,  p^»  on  the 
.subarc..  Consequently,  x  and  p  are  also  functions  of  these  quantities. 
Differentiation  of  Equations  (2.  76)  and  (2.  77)  with  respect  to  a  typical 
initial  value  a.  followed  by  an  interchange  of  the  order  of  differentiation 
on  the  left-hand  sides  (which  is  permissible)  gives 


^  d x  =  7^  jj  d x  a.  tt  .d,.P. 
dt  Sa  px11!  8  a  pnl  3a 


JL  IUlL 

dt  8  a 


-  v2h1 

x  1  da 


-  V2  hi  . 
xp  18a 


(2.  80) 


(2.  81) 
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Equations  (2.  80)  and  (2.  81)  are  a  set  of  2n  linear  first-order  homogeneous 
differential  equations  w  ith  time -varying  coefficients.  They  have  a  maximal 
set  of  2n  independent  (column)  vector  solutions.  Only  n  of  these  solutions 
are  of  interest  for  the  control  optimization  problem:  namely,  those  . repre¬ 
senting  .partial  derivatives  with  respect  to  the  initial  values  of  the  multi¬ 
pliers  pQ  These  solutions  are  represented  in  matrix  form  by 


V 

Ro 


V  p 

P0 


(2.  82) 


where  V.  x  and;  v  p  are  both  nxn  matrices.  Linear  combinations  of 

i  p  p .  -  ... 

*o 

the  solutions  (2,  82.)  of  the  form 


6x,(t)  -  v  x(t)dpQ 
-  o 


(2.  83) 


6p  (t)  *  v  p  (t)  dp 

*o 


(2.  84) 


for  small  dpQ  represent  all  possible  neighboring  extremals  abOut  the 

extreme!  whose  initial  conditions  are  x0  arid  p0  at  t  =  0  (holding 

x  fixed).  The  equations  of  these  extremals  are 
o 


x(t,  pQ  +  dpQ)  =  x(t,  pQ)  +  V  x(t)  dpQ 

Fo 


(2.  85) 


p(t,pQ  +  dpQ)  =  p(t,  pQ)  +  V  p(t)  dpQ 

1  o 


(2.  86) 
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Initial  conditions  for  the  solution  (2.  82).  of  the  differential  Equations  (2.  80) 
and  <2,  81)  must  be  determined.  For  the  first  subarc,  these  are 

%  x  (0)  =  Z.  ,  V  p  (0)  «  I  , 

pO  °0 

since  the  initial  conditions  on  the  path  are  t  =  0  ,  and 

x  (0)  =  x.  ,  p  (0)  =  p  ... 

p  o 

Differentiation  of  (2.  88)  with  respect  top  results  iri  (2;  87). 

The  initial  conditions  at  the  junction  point  t1  of  two  subarcs  are  deter¬ 
mined  from  continuity  considerations.  These  require  ;{see  Appendix  A) 


x(V  Po1 

=  x(t1  „  X ^  ,  p1 ) 

P  <t1  ,  P0> 

-  p  (t t  ,  Xj  ,  P3  );  . 

Differentiating  and  noting  that  V  x  (t.) '  =  I*  v  p  (t  )  = 

xr  1  Py 

V  p.(t1)  =  0,  and  V'  p  (t  )  =  I  then  gives 

X i  1  P-1 


7pox(ti)  +[i'<v  -  i+(tlJ 

7d  p<y  +  p"(t  -)  -  p+(tj  v  t  . 

VO  Pq 


0  , 


Here,  the  matrices  V  x(t,)  and  V  p(t,)  are  terminalvalues  for  the 

p  1  p  y  1 

.°-  _  Q  +  + 

previous  subarc,  and  x  (y  ,  p  (y  and  x  (y,  p  (y  are  limits  of 

the  time  derivatives  from  the  left  and  from  the  right,  respectively.  If  the 

control  functions  are  continuous,  the  partial  derivatives  are  continous  in 

time  at  t  ;  if  discontinuous,  the  row  vector  V  t  must  be  computed. 

1  po  1 


(2.87) 

(2.88) 

(2.  89) 

(2.90) 

(2.91) 

(2.  92) 
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It  is  determined  from  either  (A.  28)  or  (A.  35)  of  Appendix  A,  depending  on 
whether  a  constraint  is  being  added  or  subtracted.  If  (A28)  is  the  subarc 
terminal  surface,  it  is  found  from  the  implicit  function  theorem  that 


v*  ■  •  r r 


\G1  VV  + 


where  is  defined  by  the  equation  following  (A.  28) . 


The  differential  Equations  (2,8 0)  and  (2.  81)  are  usually  integrated  right 
along  with  the  set  of  Equations  (2.  76)  and  (2.  77).  In  this  manner,  the 
time-varying  coefficients  of  (2.  80)  and  (2.  81)  are  easily  calculated  at 
each  point  along  the  path. 


The  terminal  surface  of  a  path  is- usually  described  by  one,  but  not 

necessarily  all,  of  the  components  of  (2.4)  being  zero.  On  .the  surface, 

the,  first  partial  derivatives  are  evaluated  from  Equations  (2.  59)  -  (2.62). 

Examination  of  these  equations. now  shows  that  only  V  x  ,  ,m,s  still 

p  n+ 1 
0 

undetermined'.  This  maybe  evaluated  in  one  of  twoways.  The  first  is  to 
add  the  differential  equation 

d  axn  + 1  Sx  SP 

dt  da  x  o  da  p  o  da 


To  the  set  (2.  80)  and  (2,  81)  and  to'  integrate  it  along  with  the  others,. 
The  initial  conditions  for  the  first  subarc  and  at  breakpoints  are, 
respectively. 


7p  xn  +  l(0)  *  0 

0 


v  X  v 

p  n,+ 1 
*o 


-  7  X  .  (t, )  +  X  ,  .  V  t  , 

i  P0  n+l  1  n+1  P  1 
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(2.  93) 


(2.  94) 

(2.95) 

(2.96) 


The  second  .method  of  evaluating 


v  x  ,,  stems  from  a  theorem  m  the 

p  n  +  1  . 

*o 


calculus  of  variations  concerning  the  .differential  of  xn  +  j  for  a  oner- 
parameter  family  of  extremals  (Ref.  22,  p.  237  ).  Here, 


■.dXn  +  1(T>  = 


(f  +  P'f) 


o 


dT  -  p'(T)  dx(T) , 


Substitution  of 

d  x(T)  -=  x  (T)  dT  +  V  x(T)  dp 

P0  ° 

into  (2.  97),  and;  rearranging,.  then  gives 

dxn+1(T).  =  fQ  dT  -  p'(T)  V  x(T)  dp  , 

*o 

from,  which  it  follows  that 

VP  xn+  l^)  =  \  x(T)  * 

*o  ^o 

Integration  of  (2.  94)  and  evaluation  of  (2.  100)  give  the  same 
result  theoretically.  However,  for  numerical  accuracy  it  is 
better  to  perform  the  integration  than  to  evaluate  the  identity. 
Equation  (2. 100)  is  useful  for  program  check-out  and  for 
checking  the  integration  accuracy,  since  it  is  true  for  all 
values  of  the  independent  variable. 


(2.  97) 


(2.  98) 


(2.  99) 


(2.  100) 
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The  second  derivatives  come  from  application  of  the  second  variation,  see 
Ref,  22,  p.  226.  For  the  control  optimization  problem,  this  may  be  written 


T 


d2J  =  d2V  -  [ p’(T)  x(T)  dT2  +  2  p‘(T)  6x(T)  drf]  +  ;  2udT  (2. 101) 


where 


V  =  g  +  e*  \j/ 


(2.  102) 


d2V  = 


B2Y 


BT 


-,2  x 


BY’ 


B  T 


*\x(T)  +  x‘(T)  V  Vx(T) 

X 


dT 


+  2 


BY 


„2. 


V  T— r  +  X»(T)  V  V 
X  BT  X 


(2. 103) 


5  x(T)  dT  +  ox'(T)  7  Y  8x(T),  , 

x> 


and 


2w  =  quadratic  form  in  second  partials.  of  (2. 10)  in  all  variables* 

The  definition  of  the  6 -operation  is  given  on  P.  195  in  Bliss* 'book.  The 
expression  (2. 101)  assumes  that  the  transversality  conditions  (2, 19)  and 
(2.  20)  are  satisfied,  so  it  is  generally  an  approximation  during  the  iterative 
process.  The  terms  required  to  make  it  exact  are  usually  too  complex  to 
compute.  They  are  neglected  here  on  the  grounds  that  the  transversality 
conditions  are  nearly  satisfied,  and  the  Newton -Raphson  method  is  known 
to  be  insensitive  to  errors  in  the  coefficient  matrix. 


When  the  comparison  curves  are  restricted  to  the  neighboring  extremals 
(2,85)  and  (2,86),  the  integral  of  (2. 101)  can  be  reduced  to  (Ref.  22,  p.  245) 


2  w  dt  =  -  6  p»(T)  6  x(T)  , 


(2.  104) 
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where  5p  and  6x  are  defined  by  (2.  83)  and  (2,  84),  respectively. 
Substitution  into  (2. 101)  then  gives 


o  n  T  F  •• 

d  J  =  dT,  dpj}  :  W 

.  J- 


v  W 

:D- 

1  O 


(2.  105) 


,  v«  w  7*  w  dp 

p  p  o 

*o  o 


where  the  (usually  approximate)  partial  derivatives  are  comjputed  from 


W  =  -  p’(T)  x(T)  + 


52  V 


3-V  F  9 

+  27  —  lx(T)  +  x»(T)  V“V  x(T)  „  (2.106) 


3  T2  x  \  a  T  I 


sv\  ,  2  h 

V  W  =  -  p’(T)  V  x(T)  +  V  - 1+  xe(T)  7  V  V  x(T)  ,  (2.107) 

V  Pa  x  UT  x  ..  p. 


V2  W  = 


V«  p(T)  7  x(T)  +  V»  X(T)  VZV  V  x(T)  .  (2.108) 


x  P, 


sym 


Note  that  (2.  100)  may  be  substituted  into  (2.  107)  and  that  the  partials  simplify 
considerably  if  (2.  101)  is  a  linear  function. 


It  was  tacitly  assumed  above  that  breakpoints  contributed  nothing  to  either  the 
first  or  the  second  variation.  This  point  can  be  checked  by  splitting  the 
integral  of  Equation  (2.  1)  into  several  integrals,  one  for  each  subarc,  and 
expanding  each  through  second-order  terms.  It  is  found  that  breakpoints 
contribute  nothing  to  the  first  variation  because  of  continuity  of  all  the  elements 
across  such  a  point.  The  second  variation  contributions  are  zero  if  the 
breakpoint  is  not  a  corner.  Corner  points,  apparently,  require  further 
examination  for  individual  problems. 
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$ 


[ 


T 


where  T  is  the  value  of  t  on  thv-  terminal  surface;  The  Hamiltonian  for 
this  system  is 


H  =  PjXg  +  p2  ("'Xj  +  uj  +  u 
The  Euler- Lagrange  equations  are,  thus, 


•P*  =  P2 


“  P- 


P2  +  2  u  =  0  . 


The  optimal  control  function;  from  (3.9)';  is  seen  to  be 


u  =  - 


This  control  function  automatically  satisfies  the  minimum  principle -and 
the  Clebsch  necessary  condition.  Equations  (3.  7)  and:  (3/8)  are  linear 
equations  with  constant  coefficients  . and  may  be  solved  directly.  In  terms 
of  initial  conditions. 


rpx(t) 

‘ 

cos  t 

sin  t 

’pi " 

0 

p2(t> 

-  sin  t 

cos  t 

p2 

_  _ 

0 

— — 

where  p1  and  p«  are  the  initial  conditions  on  p  and  p0. 
oo  ^ 

The  homogeneous  part  of  Equations  (3. 1)  and  (3.  2)  is  the  same  as  that  for 
Equations  (3.  7)  and  (3.8).  The  system  is  thus  self-adjoint.  The  solution 


(3  .  6) 


(3.  7) 

(3.8) 

(3.9) 


(3. 10); 


(3.  11) 


! 
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of.  Equations  (3. 1),  and  (3. ,2)  can.  then  be  written  as: 


x^t) 

x2(t) 


cos  t  sin  t  x.  \ 

6 


-sin  t  cos  t  :  x 


(3S 12) 


cos  t  sin  tl  c  r  cos  t  -sin  F  0 


-sin  t  cos  -y  J  j^sin  t  cos  tj  j^l  J’ 
-  o  : 


P2  (t)  d,T 


After  substitution  and  integration, 


Xj(t) 


X2W 


cos  t  sin  t  x 


:sin  t  cost  lx 


(3.  13) 


(t  cos  t  -sin  t)  t  sin  t 


-t  sin  t  (sin  t  +  t  cos  t) ,  p 


Finally,  the  evaluation  of  expression  (3.  3)  is 


12  2 
J(t)  =  -g-  (t  -  sin  t  cos  t)  -  2  p^  p2  sin  t 

O  0  0 


+  p0  (t  +  sin  t  cos  t) 


(3.  14) 
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On  the  terminal  surface,,  then. 


1  _  .2; 


J(T)  =  —  (T  -  sin  T  cos  T)  r  2  px  p2  sin“  T 


o  o 


*  pr.2  (T  +  sin  T  cos  T).  * 


'  xt(T) 


ix2(T) 


cos  T  sin  T  x 


'-.sin  T  cos  T  ,  x2 


(T  cos  T  -  sin  T) 


•  ■  v 

>  o 


-  T  sin  T 


T  sin  T 


(sin  T  +  T  cos  T)  p2  ; 

: ;  o  a 


(3,  16)< 


The  function  F  of  Equation  (2.  28)  for  the  equivalent  minimization  problem 
may  be  written,  from  (3.  15),  (3.4)  arid  (3.5),  as: 


F  =  J(T)  +  e1(T  -  K)  +  e2  x^  (T) 


(3.  17) 


The  first  necessary  condition  then  gives 


3 F  1  o 

— ~~  =  o  =  —  (p,  sin  T  -  p0  cos  T)  +  e.  +  e0  /-  x  sin  T 
ST  4  o  o  1  M 


(3.18) 


+  x  cos  T  +  —  p.  T  sin  T  -  p0  (sin  T  +  T  cos  T) 


4  o 


BF  1 

spi"  p  '  T 

O 


2 

p  (T  -  sin  T  cos  T)  -  p„  sin  T 

O  "  o 


e  • 


(T  cos  T  -  sin  T) 


3  F 
BPo 


*  0  = 


1 

T 


i 


-  p.  sin^  T  +  p2  (T  +  sin  T  cos  T) 
Ao  o 


(3.19) 


(3.20) 


”  4 


T  sin.  T  ... 


These,  together  with  (3.  4)  and  (3.  5)  are  to  be  solved  for  T,  p^ 

< 

P2  >  ej  and  e2  *  The  solutions  are 
o 


K 

A  cos  K 


(3.  21) 
(3.22) 


A  sin  K 


(3.  23) 


Aj"x^  K  sin  K  -  x2  (K  cos  K  -  sin  K) 
L  o _ o _  J 

(K  -  sin  K  cos  K) 


(3.  24) 


A  = 


jx^  cos  K  + 
L  o 


x2  sin  K 
o 


] 


(K  -  sin  K  cos  K) 


(3.25) 


provided  K  ^  0  . 


40 


,jvr*?ti*** 


The  point  defiried  b^  Equations  (3.  21)  through  (3  .  25)  is  truly  a  minimizing 
point.  This  will  now  be  shown  by  constructing  the  matrix  of  Equation  (2.  36), 
The  set  (3.  4)  and  (3.  5)  is  first  completed  by  adding,  the  equation 


^3  =  V®  +  y‘ 

Then  the  Equation  (2.  3  5)  may  be  written,  for  the  critical  point,  as 


r 

i 

0 

0 

< 

I  AT 

>  -- 

■  0 

1  *1 

*11 

*12' 

Api  - 

■  =  ■ 

0 

<  j 

0‘ 

, 

*2 

X21 

*22 

' 

,  *P2 

Ay 

V 

0 

/ 

— .  - 

where 


x 


11 


dx 


1 


o 


12 


aX2  x 
SPj  '  22 


Equation  (3.  27)  may  be  rewritten 


7  AT 

■  1  0  0 

"c 

Ap  f 

-  = 

a21  a22  a23 

0 

0 

- 

Ap2 

a31  a32  a33 

’  Ay.. 

0 

—  mi 

_ 

where 


a21 

C(x12x2 

*22*2  ^  a22 

Cx22' 

a23  =  “  Cx: 

a31 

C(x21il 

*11*2^'  a32 

C*21' 

a33  =  Cxll 

1 

~C  = 

(xl-l  X22 

*21  *12^  ' 

<3.  26) 


(3.  27} 


(3.  28) 
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Also,  at  the  critical  point. 


V  2f 


0 

4  K  sin  K 

4 

.(sin  K  +  K  cos  K) 


4  K  sin  K 

4 

(K  r  3in  K  cos  K) 

•  2  __ 
sin  K 

'  4*  ‘  ‘ 


-  4<sin  K+Kcos  K)^ 
sin2  K 

4 

t 

~(K  +  sin  K  cos  K) 


(3.  29) 


After  all. the  substitutions  have  been  made,  it  is  found  that  the  matrix  A^, 
of  Equation  (2.  36)  is  the  scalar 


_  (K  -  sin  K  cos  K) 
4  (K2  -  sin2  K  ) 


Since  K  >  sin  K  for  all  K  >  0,  A^  >  0.  It  is,  thus,  concluded-  that  the 
solution  is  indeed  a  minimizing  solution,  since  there  is  only  one  critical  point. 


The  accessory  differential  Equations  (2.  80)  and  (2.  81)  for  this  problem  can  be 
written 


d 

W 


Sxl 

TsT 

3x2 

B  a 

BP1 

Ba 

5Pg 

d  S- 

_  — 

-1 


0 
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0 

0 

0 


0 

0 

0 

-1 


0 

j— - 

I*5 

1 

8x2 

"2 

B  a 

1 

B  P| 
Ba 

0 

its 

Sa 

(3,  30) 
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It  can  readily  be  verified  that  the  solution  of  C3,  30)  is 


Bx1  dx^ 

■ 

sint  -  t  cos  t  -  t  sin  t 

'  ' 

o 

sp2 

0 

; 

4  4 

;  ax2 

ox2 

, 

t  sin  t  -  (sin  t  +  t  cos  t) 

:  5pi 

0 

b,p2 

O' 

.  4  4 

7  P= 


Spl> 

-pi 

opr 

sp2 

Q 

0 

bp2 

Sp2 

Sp, 

3P.2. 

cos  t  sin  t 


-  sin  t 


cos  t 


,  V  p  (0)  =  I 
-  o 


The  first  derivative  of  (3.3)  with  respect  to  pQ,  from  (2.  100),  is 


V  J  =  -  pMT)  V  x  (T)  =  ]. 

P  P  4 

-o  *o 


2 


;(p  (T  -  sin  T  cos  T)  -  p2  sin  T) 
Lo  o 


.  2 


(-  p1  sin  T  t  p2  (T  +  sin  T  cos  X)) 


x(0)  -  Z 


(3.  31) 


(3.  32) 


(3.33) 


It  is  seen  that  this  corresponds  to  the  results  of  differentiating  Equation  (3.  14). 


Since  the  constraint  Equations  (3.  4)  and  (3*  5)  are  linear,  it  is  seen  that  the 
matrix  (2...  105)  of  second  partials  reduces  to 


V-  2  J  *  V »  x(T)  v  p(T)  = 


Q 


o 


T  r  sin  T  cos  T 


■  sin-2  T 


.  2 

-  sm  T 


T  +  sin  T  cos  T 


4 


<3,34) 


This  agrees  with  the  corresponding  terms  of  (3,  29).  when  the  identification 
T  -  K  is  made. 


THE  HIGH  LIFT-DRAG  RATIO  VEHICLE  OPTIMIZATION  STUDY  (B) 

The  Problem  Statem  ent  ,  (B  - 1 ) 

A  trajectory  is  sought  which  minimizes  the  functional  (4,  6)  for  a  fiat-plate 
vehicle  with  aerodynamic  coefficients  represented  by  Equations  (1.  3).  The 
data  for  these  formulas  is  C^q  =  1.  82,  C^q  =  0.  042,  and:  =  1.  40, 

The  motion  is  assumed  to  .be  governed  by  equations  of  the  form  <1, 1)  and  the 
atmosphere  is  specified  by  the  relation 


P  =  P0  e 


-0  R§ 


(3,  35) 


The  terminal  surface  is  represented  by  the  single  equation 


-•  T  -  K  =  0  , 


(3,36) 


where  K  is  a  specified  constant.  The  pilot5 s  acceleration,  computed  from 
the  formula 


=  it  P-Y.1-  „  lei  +  C? 


P  2  m  g. 


D 


(3,37) 
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Xt- 

is  constrained  by  the  inequality 

B  -  a  =  0  ,  (3.38) 

fP 

{The  constant  B  is  generally  taken  as  lO  .g. )  Initial  conditions  for  the 
path  are 

v  -  35,  000  feet,  per  second 

o  ■ '  ■ 

7Q  ~  ~  5  degrees 

h  =  400,  000  feet 

o  1 

C  -  o  . 

°o 

The  criterion  (Is  6)  is  chosen  to1  correspond  to  a  vehicle  which  radiates 

heat  most  efficiently  at  a  specified  temperature.  The  constant 

# 

Q  =  195  BTU/sq.  ft. /second  fixes  this  temperature  around  30Q0°R 

for  this  particular  vehicle,  and,  since  the  trajectory  is  selected  by  the 
criterion,  the  heating  rate  of  the  vehicle  deviates  in  the  mean- square 
sense  as  little  as  possible  from  this  constant.  To  simplify  the  computations, 
the  heating  law  is  taken,  following.  Chapman  (Kef.  8)  as 

q  =  Cpl/2  w  ,  (C  =  2  x  10'8)  ,  (3.39) 

A  more  realistic  formula,  which  includes  a  radiative  heating  term  in 
addition  to  this  formula  representing  cohvective  heating,  is  used  in  the 
next  subsection. 

Since  T  is  a  constant  by  (3.  36),  it  may  be  omitted  from  the  criterion; 
so  the  integrand  f  of  Equation  (2. 1)  may  be  identified  with  the  integrand 
of  (1.8): 

f0  =  (Q  -  q>2  .  (3.40) 
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The  function  g  of  (2. 1}  is  taken  as  zero.  Thus,.  Equation  (2, 10)  for  this 
application  reads 

F  =  f  +  p  (f  -  x)  +  p.  (B  -  a  -  a2 )  ,  (3. 41) 

o  p 

where  the  terms  on  the  right-hand  side  correspond  to  relations  {3.40), 

(i.  1),  and  (3.  38),  respectively.  The  range  coordinate  Q.  does  not 
appear  explicitly  in  these  terms  and,  hence,  is  a  cyclic  (or  ignorable) 

coordinate.  The  corresponding  Euler-Lagrange  Equation  (2.1,1)  is.  i-  « 

• 

p^.  «  0;,  hence,  p4  must  be -a  constant}  This  constant  is  zero  in 
order  that  the  transversality  condition  (2.  20)  is.  satisfied  on  the  terminal 
surface  (3,  36).  Hence,  the  fourth  equation  of  system  (1.  1)  may  be 
eliminated.  Further  evaluation  of  the  transversality  condition  gives 

PX(T)  =  p2(T)  =  p3(T)  =  0  , 

and 


H  =  e  , 

where  e  is  a  constant  which  is  yet  to  be  determined.  To  compress  the 
notation,  the  remaining  three  coordinates  are  written  as 

Xj  -  v  ,  x2  =  y  ,  x3  =  ^  . 


The  Unconstrained  Subarc  (B-2) 


The  multiplier  p  in  expression  (3.41)  is  zero  on  any  subarc  which  is  not 
limited  by  equality  in  the  constraint  relation.  This  follows  from  (2. 13) 
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t 


This  leads  to  the  equation 


tan  a 


-  ai 


$ 


where 


3Cdl  fll 

2Glo  H  ■ 


Note  that  the  sign  of  tan  a  is  determined  by  the  choice  of  the  ±  sign  in 
Equation  (3. 46) 


The  lift-drag  polar  is  traversed  once  as  a  ranges  through  tt  radians. 
It  is  thus  advisable  to  limit  the  control  function  to  a  range  of  TT  radians 
to  avoid  double  values.  Furthermore,  the  points  a  '=  0,  mare 

^2  p 

singular  points,  since  - — s-  (the  determinant  Rr  of  (A.  13)  in 
Appendix  A)  is  zero  at  these  points.  The  range  0  <  a  <  n  is 
chosen  here,  since  then,  !  sin  a  I  =  sin  a  ,  and  the  proper  sign  for 
the  equation  is  most  easily  chosen.  The  singular  points  are  removed 
by  defining  a  =  0  at  these  points.  Thus,  the  range  of  a  is 

0  £  a  <  tt  . 

The  proper  choice  of  sign  is  determined  by  considering  the  minimum 
principle  which  yields  the  inequality 


a  •- 


2 

sin  a  |  •-  Cql  PjV  sin  a 


2 

<:  sin  A  -  p^v  sin  A 


+  CLOP2sino 
*  CLO  p2  sin  A 


cos 


a 


cos  A 


]. 


where  A  is  any  value  of  the  control  function  satisfying  Inequalities  (3.  47) 
and  (3.38).  Substitution  of  (3.45)  into  (3.48)  and  rearranging,  then  gives 

Pg  sin  a  tana  £  pg  sin  A  tan  A  , 
providing  tan  A  satisfies  Equation  (3.  45)  . 


(3,46) 


(3.47) 


(3.48) 


(3.49) 
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Since  sin  or  *  0  arid  the  sign  of  tan  or  is  determined  by  the  ±  sign  of 
Equation  (3.  46),  it  is  seen  that  (3. 49)  is  satisfied  by  the  choices: 


If 

if 


p9  >  0,  choose  -  sign  (-S-  <  a  <  tt  )  ; 


Pg  <  0,  choose  +  sign  (  0  <  a  <  . 


The  behavior  where  ~  0  is  examined  by  a  limiting  process,  since 
Pi*  ^2  anc*  v  are  con^nuous  functions  of  time.  The  binomial  series 
expansion,,  assuming  large  a,  gives 


t  2 


=  I  a  !  + 


If  a  >  0-  and  the  minus  sign  were  chosen,  then 

tan  a  =  lim  -  2  |  a  |  ,  d  -*  (from  the  left)  . 

a  -*« 

If,  a  <  0  and  the  plus  sign  chosen, 

tan  a?  =  lim  2  la  I  ,  a  —•  (from  the  right)  . 

a  -»  co 

On  the  other  hand,  if  a  <  0  and  the  minus  sign  chosen, 

tan  a  =  lim  -  — =  0  ,  d  -♦  rr  (from  the  right)  , 

a  -*  «  }a  1. 

Finally,  with  a  >  0  and  the  plus  sign, 

tan  u  =  lim  =0,  d  -*  0  (from  the  left)  . 

a  -*«  1  a  I 
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The  list  two  cases  show  that  it  is  possible  for  at  tobediscontinuous.  For 
example,  if  p2  passes  through  zero  and  at  was  near  tt  ,  then  at  will  jump 
to  zero  as  p  .  ;goes  through  zero.  This  difficulty  is  cleared  up  if  the  proper 

A  ' 

interpretation  is  made.  Since  tan  (tt  +  a)  =  tana,  a  may  be  thought  of 
as  being  continupus  across  jumps,  but  in: the  dis-allowed5 region  tt  js  a  <  2tt. 
The  control  function  discontinuity  will  be  retained  since  computer  results  can 
properly  be  interpreted,  and  since  it  allows  I  sin  a  !  -  sin  a  in  all  equations. 
Finally,  although  of  little  practical  importance,  if'  and  p2>  are  simultaneously 
zero  at  a  point,  the  limit  ratio.  p^/p2  is  Pj/p2  at  thatpointby  L'Hospital’s 
rule. 

To  sum  up,  the  reduced  differential  equations  of  the  extremals  for  subarcs 

haying  a^  <  B  are  the  first  three  of  (1. 1).  and  system  (3.  42)  .  The  control 

function  is  computed  from  Equation  (3.  46)  in  the  range  0  £  a  <z  tt  according 

to  the  rules  following  Equation  (3. ,49)  .  The  subarc  terminal  surface  is 

a  =  B,  provided  that  a  5^  0  . 

P  P 


0  = 


af,  bf0 

1  t  Pr  1 


i,a  a 


'2  a  OL 


-  H 


a  aF 


(3.  51) 


The  control  function  is  determined  from  a  =  B  and  the  multiplier  p 
from  (3.  51).  Expansion  and  rearrangement  of  the  constraint  equation 
B  gives 

H  ' 

v2  - 


_2  .  3 

CDO  '  Sln 


Of 


.  3 

cx  sin  « 


+  Cg  Sin  a  + 


(3.52) 


where 


b 


2m  BgQ 

-r—  -j- 

S  p  v 


o  =  c  -  c 

cr  -  DL  VLO 

C2  =  CLO 

C3  =  2  CDO  CDL  * 

The  right-hand  side  of  (3.  52)  is  a  function  only  Of  the  vehicle  aerodynamic 
coefficients  (constant^)  and  sin  a  .  It  is  zero  when  or  =  0  or  tt  ,  and 
maximum  when  a  =  -y  .  From  this  it  follows  that  b  must  satisfy 


CDO  *  b  s  *CDL  +  '*  *3*  53  ^ 

The  physical  interpretation  is  that  if  >  b  at  any  point  of  the  path, 

then  < a  has.  gone  to  zero  in  the  futile  attempt  to  keep  a  =  B.  The  inte- 

P 

gration  is  stopped  if  condition  (3.  53)  is  not  satisfied. 
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Newton* s  method  was  chosen  to  extract  a  from  (3.  52)  .  The  iteration 
equation  is 


A  a  "= 


/.  2  ^,2  v 

(b  -  C  -  sin 


sin  a  cos  a  6  c.  sin 

1 


.  * 

3  3  ! 

i  a  : ’  sin  a  +  Cg  sin  o'  +  Cg 

» I  -  - .  -  -  "»' 

sin  a  +  4Cg  sin  or  +  3  Cg  • 

m  ' 


(3.5# 


A  Aa  is  calculated  from  (3.  54),  using  an  assumed  a  (usually  the  last 
integration  step  value)  This  is  added  to  a  ,  and  the  hew  or  is  used  to 
calculate  a  new  A  a  ,  The  process  continues  until  A  or  is  negligible. 


In  the  vicinity  of  or  =  0,  j  or  n  ,  the  denominator  of  (3.  54)  is  likely 
to  be  quite  small.  To  avoid  this  difficulty,.  Equation  (3.  52)  was  expanded 
as  a  series.  The  results  are: 

a  3  ^*.c*0\1,, 

If  a  ~  0  ,  (sin  or  =  a  -  -—) ,  then  a  1  ; 


3 

ot  ~  tt  ,  (sin  or  =  (tt  -  a)  -  ^  ^  )  , 


then  a  =  tt  - 


b2  -  C  2  \  1/3 
^DO 


a  ~  ,  (cos  a  =  (-^  -  a  )  -  -  a  )3  ) , 


then  1“  9  9  H  1/2 

n  ±  <CDT.  *  Cno'  -  ^ 

3  CDL  (CDL  +  CDQ)  -  CLQ 


where 


the  +  sign  is  used  if  Pg  >  0  ,  and 

the  -  sign  is  used  if  pg  <  0  . 
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If  Pg  =  0  ,  choose  sign  so  that  tan  a  is  continuous.  The  choice  of 
signs  again  comes  from  the  minimum  principle. 

-I 

It  is  also  noted  that  the  determinant  R9  of  Equation  (A.  23)  in  Appendix  A 


Since  this  must  never  be  zero  when  a  =  B,  a  must  never  be  6,  r~- 

p  ■'  •  2 

of  TT:  . 

Finally,  jjl  is  calculated  from  Equation  (3.  51)  in  the  form 


The  reduced  differential  equations  of  the  extremals  for  the  constrained 
subarc  are  now  the  first  three  equations  of  (1.  1)  and  the  system  (3.  50). 
The  control  function  (0  <  a  <  -Tp  or  <  a  ^  <  rr)  from  either  the 
iterative  Equation  (3.  54)  or  the  small  angle  equations,  and  the  multiplier 
H  (s:  0)  comes  from  (3.  55)  . 


The  Subarc  Junction  Points  (3-4) 

It  is  known  that  the  multipliers  p  and  the  Hamiltonian  H  are  continuous 
at  the  junction  point  of  constrained  and  unconstrained  subarcs.  It  is 
necessary  for  this  problem  that  the  points  at  which  a  =  0  ,  -y  be 

ruled  out  as  junction  points,  since  the  determinants  and  of 
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Equations  (A.  13)  and  (A.  23)  must  be  non- singular  at  such  points.  For  ail  other 
values  of  or  ,  it  turns  out  that  the  control  function  must  be  continuous  at 
junction  points.  This  is  readily  verified  from  examination  of  the  Hamiltonian. 

It  is  an  analytic  function  of  x,  p  and  or  .  Since  x  and  p  are  continuous  at 
junction  points,  and  since  the  Jacobian  Rj  is  non- singular,  it  follows  that 
or  is  continuous  at  such  points.  Similarly,  since  Rg  is  non- singular. 

Equation  (3.  55)  shows  that  \i  must  be  continuous  at  such  points.  Thus,  > 
must  start  and  end  with  the  value  zero  on  constrained  subarcs.  It  then 
follows  that  the  constrained  .subarc  terminal  surface  is  p  -  0,  provided 

(x  ^  0  at  the  junction  pointi 


The  Newton-Raphson  Equations  (B-.5) 


The  form  of  the  equivalent  minimization  problem 
Minimize  the  function  J(T,  pQ);  corresponding  to 
the  variables  y*  =  {T,  p1Q,  p2Q,  ,p3Q)  subject  to 


considered  here  is: 
the  functional  (1.  6)  in 
the  constraint  equation 


^  =  T  -  K  =  0  . 

The  function  F  of  Equation  (2.  28)  is 


(3.  56) 


F  =  J(T,  Pq)  +  e(T  -  K)  ,  (3.  57) 

which,  when  differentiated  with  respect  to  T,  gives 

“tf^  P0‘  e)  0  =  j  +  e  +  dpQ  ,  (3#  5g) 


since  the  constraint  (3.  56)  implies  dT  =  0,  Tie  constant  e  of  (3.  58)  can 
always  be  chosen  to  satisfy  the  identity.  Differentiation  of  (3.  57)  with 
respect  to  the  other  variables,  setting  the  result  to  zero,  and  a  little  re¬ 
arrangement  then  gives  the  Newton-Raphson  iterative  equation 

dp  =  -  (vj?  J)*"1  V.  J  .  (3.  59) 

°  po  po 
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>x_.' _  •  <nms .  -  .  , 


2 

The  sufficiency  condition  reduces  to  showing  that  the  matrix.  J  is 

'*  Q 

positive -definite.  This  follows  since  in  the  quadratic  form 


dy*  VyF  dy  >  0 


(3.  60) 


dy.  s  0  and  the  other  components  of  dy  are  arbitrary. 


The  steepest  descent  equation , 


dp 


d  J 


o 


v  J  j 


v»  J 
?o 


(3.  61) 


is  readily  verified  using  the  methods  of  subsection  11(D).  The  magnitude 
of  dJ  is  chosen  by  experience,  usually  smaller  and  smaller  as  the  optimum 
is  approached. 


Computer  Results  (B-5) 

Some  of  the  computer  results  are  displayed  in  Figures  3-1  through  3-4. 
These  trajectories  were  obtained  using  the  method  of  estimated  partial 
derivatives  of  subsection  II  (E)  . 

The  terminal  times  for  the  optimal  trajectories  in  Figures  3-1  to  3-  3  are 
200,  440  and  550  seconds,  respectively.  It  is  noted  that  the  angle  of  attack 
histories  are  roughly  the  same,  going  from  small  to  large  values  after 
about  100  seconds  of  flight.  This  corresponds  to  a  shift  from  maximum 
positive  to  maximum  negative  lift.  It  is  further  noted  that  this  shift  takes 
place  near  the  bottom  of  the  first  pull-out,  where  the  heating  rate  is  at  its 
peak.  The  vehicle  apparently  does  not  possess  enough  lift  to  avoid  the 
heating-rate  peak  for  the  given  initial  conditions,  nor  can  it  avoid -the 
mild  skip  shown  in  the  figures. 


Figures  3-2  and  3-3  show  that  for  larger  termir.V  tiroes  the  vehicle  dtve.& 
rather  sharply  hear  its  terminal  points.  Evidently  the  integral  of 
Equation  (1.  61  is  minimized  by  this  dive,  in  spite  of  the  large. rise,  of 
short  duration,  in  the  heating  rate.  It  should  be  remembered  that  no 
restrictions,  other  than  fixed  terminal  time,  were  placed  on  the  terminal 
point. 


Figure  3-4  is  included  to  illustrate  the  sensed  acceleration  constraint.  It 
does  not  represent  an  optimal  path.  The  angle-of-attack  history  is  similar 
to  the  other  results  until  the  constraint,  set  at  10  g’s,  becomes  an  equality 
constraint.  From  this  point  onward,  the  angle  of  attack  changes  such  that 
10  g's  is  never  exceeded.  The  magnitude  of  the  additional  multiplier  p 
is  shown  at  the  bottom  of  the  figure. 

Do  the  results  of  Figures  3-1  to  3-3  represent  relative  Or  absolute  .minimum s 
for  the  given  problem  ?  An  attempt  to  answer  this  question  was  made  by 
re-optimizing  with  very  different  initial  conditions  for  the  multipliers. 
Computations  were  carried  far  enough  to  show  that  the  same  optimal  path 
would  have  resulted,  indicating  that  the  paths  probably  represent  absolute 
minimums.. 

The  importance  of  these  results  is  not  that  they  are  optimal  paths,  for  they 
have  ;many  practical  short-comings,  but  that  they  confirm  the  validity  of  the 
automatic  optimization  scheme.  Further  improvement  of  the  optimization 
method  was  judged  to  be  more  significant  and  important  at  this  stage  than 
the  computation  of  an  operationally  desirable  re-entry  path. 
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The  next  stage  of  development  was  the  method  for  computing  the  partial 
derivatives  explicitly  (see  Subsection  11(F)).  The  computer  program  incor¬ 
porating  this  method  calculated  the  first, par tials  directly  from  Equation  (2.  94)* 
whereas  the  second  partial  derivative  matrix,  from  Equation  (2. 108)>  took  the 
form 


Attempts  were  made  to  duplicate  the  results  of  Figures  3. 1  through  3.  3  to  see 
if  further  improvement  in  the  trajectories  could  be  made.  Some  reduction  in 
the  optimal  criteria  was  obtained,  ’  it,  in  all  cases,  the  second  derivative 
matrix  was  so  close  to  singularity  that  either  bad  prediction  of  multiplier 
changes  resulted,  or  the  matrix  became  indefinite.  This  behavior  indicated: 
that  the  surface  was  extremely  flat.  It  was  decided  to  abandon  this  program  in 
favor  of  the  more  realistic,  more  complex  and,  hopefully,  less  sensitive 
problem  considered  in  the  next  section. 

This  program  did  show  that  the  method  of  computing  partial  derivatives  was 
superior  to  that  of  estimation. '  It  also  ran  from  three-  to  five -times  faster 
than  the  previous  program  on  the  computer. 


THE  LOW*' LIFT-DR  AG  R.iTIO  VEHICLE  OPTIMIZATION  PROBLEM  (C) 

The  Problem  Statement  and  the  Euler -Lagrange  Equations  (C-l) 

A  re-entry  path  which  minimizes  the  total  stagnation  point  heating  (1;  7)  for  the 
blunt-nose  body  for  which  the  aerodynamic  coefficients  are  given  by  Equa¬ 
tions  (1.  4)  is  to  be  found.  The  vehicle  aerodynamic  constants  are  G  =  0.  88, 
Cdl  =  0.  52  and  C^q  =  -0.  505  ,  and  the  model  equatipns  are  again  the  set 
(1. 1)  with  the  exponential  atmosphere  (3.  35).  The  terminal  surface  equations  are 


\f/1  =  v(T)  -  V  =  0 

\f/2  =  h(T)  r  A.  =  0  (3.63) 

^3  =  a  m  -  r  =  o , 

with  the  constants  V  =  1650  feet  per  second,  A  =  75,  530  feet,  and 
R  =  979  statute  miles.  Thus  only  the  final  flight  path  angle  and;  terminal 
time  are  left  unspecified;  Initial  conditions  are  taken  as  vq  =  35, 000  feet 
per  second,  =  -5. 75  degrees,  'h  =  400,  OOO  feet,  CQ  =  zero.  Inequality 
constraints  imposed  are  relation  (3,  38)  and  a  bound  on  the  control  function  u, 
given  by 

u|  -  u^  >  0  ,  (3. 64) 

with  Uj  a  constant  (16  degrees  at  present).  Equation  (3.64)  was  found 
necessary  to  produce  initial  trajectories  which  neither  skipped  out  of  the 
atmosphere  nor  dived  in  too  deeply.  It  is  to  be.  sequentially  relaxed  as 
the  proper  region  of  pQ  space  is  located  during  the  optimization  process. 

The  integrand  of  criterion  (1.  7)  is  taken  as  the  sum  of  convective  and 
radiative  heating  rates : 

q  =  qc  +  qr  „  (3. 65) 


wherevthe  convective  component  is  given  by  Equation  (3,  39)  and  the 
radiative  component  by 


(3.66) 


in  which  N  =  4  feet,  the  vehicle  frontal  nose  radius.  Equation  (3.  65) 
may  again  be  identified  with  the  integrand  f  of  Equation  (2. 1),  and  the 
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function  g  is  omitted.  Thus  Equation  (2. 10)  becomes 

F  =  q  +  P  •  (f  -  x)  +  PjCu2  -  u2  r  <j| )  +  m-2  (B  -  ap  -  a2 )  , 

where  the  state  coordinate  x  now  has  four  components  with  x^  =  ?  . 
The  Euler -Lagrange  equations  are 


-Pi  = 

3q. 

3X1 

+  p, 

l£i 

“i 

+  p2 

111 

axt 

+  p3 

Jf3 

axx 

+  p4 

«4 

BaP 
^2  ax^ 

• 

5f2  * 

af3  . 

«4 

”  p2  = 

pl 

SX2 

p2 

ax2 

p3 

ax2 

p4 

3x2 

dq 

af2 

Sf4 

Bap 

.p3  = 

dx„ 

+  P1 

ax. 

+  p2 

axQ 

'  +  P4 

Sx, 

^2 

dx. 

^  P4  =  0  ,  <p4  -  p4o  ) 


I 


0  =  p 


1  du 


3f, 


+  P 


2  du 


-  2Pj u  - 


u. 


du 


The  Unconstrained  Subarc  (C-2) 


Both  the  multipliers  and  are  zero  here,  as  noted  in  Subsection  11(B). 
Then,  when  the  substitutions  from  Equations  (1. 1)  and  (1.4)  have  been  made. 
Equation  (3. 68)  becomes 


tan  u 


"CLOp2 

CDLP1V 
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(3.,67) 


(3.68) 


(3.69). 


T  f 


where  u  is  centered  about  zero  by  the  constraint  (3. 64),  i.e., 

-Uj  i  u  *  Uj  .  (3.  70) 

The  minimum  principle  equation  is 

-  Pi  vCDL  cps  u  +  p2  CLG<  sin  *  *  ‘Pl  ^DL  oos  U  +  **2'CL£>  sin  °> 

(3.  71) 


where  U  is  iuiy  admissible  value  in  the  range  (3.  70). 

The  left  hand  side  of  (3.  71)  may  be  considered  as  a  dot  product  and  the  choice 
of  a  unit  vector  (cos  u,  sin  u)  which  has  minimum  dot  product  with  the  vector 

<-P*  vCDL'  P2CLO)  is 


COS  u  = 


cdl  »i  v 


V<CIX>P2)2  *  <CDLP1V)2 


-C 


(3.  72) 


sin  u  = 


LO  p2 


V  ^CLO  P2  *  *  *CDLP1 


which  is  parallel  but  in  the  opposite  direction. 

Then  from  the  signs  of  p1  and  p2,  assuming  CL0  negative,  it  follows  that: 


If 


p2  = 
P2  > 
p2  > 
p2  > 
p2  = 
P2  < 
P2  < 
p2  < 


0 

0 

0 

0 

0 

0 

0 

0 


and  p1  > 
pi  > 


< 

< 

> 


0 

0 

0 

0 

0 

0 

0 

0 


then  u  =  0 


0  <  u  <  f 

u  =  f 

i  <  u  <  n 


u  =  ±rr  (bang  condition  if  u.  =  tt) 


(3.73) 


-  ~  <  u  <  0 

-  -  -t 

-n  <  u  <  | 


1 
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There  are  no  singular  points  if  arid  pg  are  never  simultaneously  zero. 

This  is  verified  by  computing  the  determinant  of  Equation  (A.  13).  The 
subarc  ends  either  with  (3. 38)  or  (3.  64)  becomes  zero. 

The  Constrained  Subarc  u  =  ±u^  (C-3) 

Let  cp  be  the  angle  defined  by  Equations  (3.  72)  and  the  sign  conventions 
given  by  (3.  73).  Then  substitution  into  the  minimum  principle  Equation 
(3.  71)  /»ives 

cos  (cp  -  u)  cos  ( cp  -  U)  ,  (3.  74) 

which  is  satisfied  if  u  and  cp  <  ±tt  ,  have  the  same  sign.  The  condition 
cp  -  tt  indicates  a  bang.  Furthermore,  substitution  into  (3.  68)  (with 

=  P)»  8ives 

S  p  v  rrs - - — 9  — — — 1 — -o  sin  ( cp  -  u) 

*1  "  "  1m"  V(qLQP2)  +  (CDLPlv)  — JT 

Since  u  and  sin  (cp  -  u)  have  the  same  Sign,  Uj  =  0,.  as  required  by 
Equation  (2, 14). 

There  are  no  singular  points;  this  is  easily  shown  by  computing  the  determinant 
Rg  of  Equation  (A.  23).  It  is  also  easily  shown  that  the  control  function  is 
continuous  at  the  junction  between  constrained  and  unconstrained  subarcs. 

Thus  Py  from  (3.  75),.  must  start  and  end  with  value  zero,  since  at  such 
points  u  =  cp.  Then  the  terminal  surface  is  =  0,  provided  4  0. 


(3.75). 


The  Constrained  Subarc  a  =  B  (C-4) 

.  ... . p  ~ , 

The  control  function  is  determined  from  a  =  B  and  the  multiplier 
from  Equation  (3.  55),  Expansion  and  rearrangement  of  the  constraint 
gives 

2 

a  cos  u  +  2  b  cos  u  +  c  =  0 


where 


a  =  C 


DL 


-  C 


LO 


b  CWt 

DO  DL 


c  =  c  2  +  C  ^ 
c  '“"DO  -  LO 


2  m  gQB 


S  p  v' 


It  follows  by  substituting  u  =  0,  ±rr  into  (3.  76)  that  the  inequality 


(CDO  *  CDL)  * 


2  m  g  B 


S  p  v 


2  *  <CDO  ~  CDL^  >  0 


must  hold.  Furthermore^  the  determinant  Rg  of  Equation  (A,  23)  is 
singular  at  these  points. 

Solution  of  (3,  76)  for  u  gives 

b 

cos  u  = — 
a 

where  the  omitted  root  falls  outside  the  range  I  cos  u|  <:  1.  It  is  easily 
shown  that  the  term  under  the  square  root  in  (3.  78)  is  positive  by  sub¬ 
stituting  the  upper  limit  of  (3.  77)  into  the  expression  and  evaluating. 


<3.  76) 


(3.77) 


(3.  78) 
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The  minimum  principle  again  takes  the  form  <3.  74).  except  that  this  time 

rU  <  cp  <  u.  (3.79) 

Otherwise,  cp  and  u  would  be  identical.  Then,  once  again,  cp  4  0  and 
u  must  have  the  same  sign,  and  cp  =  0  is  the  bang  condition.  It  can 
further  be  shown  that  the  multiplier  n2  becomes,  when  all  the  substitutions 
have  been  made. 


Since  sin  (u  -  cp)  and  sin  u  have  the  same  signs,  \i^  s  0  as  required. 
Again,  u  is  continuous  at  junction  points,  so  n2  must  start  and  end  at 
zero;  and  =  0  with  n2  ^  0*  describes  the  subarc  terminal  surface. 


The  Modified  Newton-Raphson  Equations  (C-5) 


The  terminal  conditions  for  this  problem  are  of  the  form  (2.  47)  which  means 

the  iterative  Equation  (2.  54)  applies  with  p2(T)  =  0,  the  single  multiplier 

constraint.  Let  "H..  and  C..,  i,  j  =  1,  ...  4,  be  the  elements  of 

V  x(T)  and  V  p(T),  respectively.  Then  the  matrix  (2.  55)  may  be  written 
?o  Pq 


v(T) 

”12 

”13 

?  (T) 

^31 

,T132 

^34 

C  (T) 

”41 

”42 

”43 

^44 

0 

fj(o) 

f2(°) 

f3(°) 

u{0) 

p2m 

^21 

^22 

^23 

C24 

(3.81) 
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The  last  vector  of  Equation  (2;  54)  assumes  the  form 


0 

:  .  C(y0) 

H(y0) 

P? 

where  the  first  component  v(y  )  -  V  =  6  is  the  stopping  condition  on  the 
integrations  and,  hence,  is  satisfied  by  every  trajectory. 

The  matrix  of  second  partial  derivatives  in  Equation  (2. 105)  is  the  five-by 
five  matrix 


v2  J  (y  ) 

y  Jo' 


-p'(T)  x  (T) 

-  V*  x,(.T)  p'  (T) 
po 


-p'(T)  V  x(T) 
Po 


V»  x  (T)  V  p  (T) 
_  po  po 


sym. 


(3. 82) 


Both  methods  of  Subsection  11(C)  were  used  to  reduce  this  to  a  two-by-two 
matrix  to  be  tested  for  positive-definiteness.  Eor  the  determination  using 
inequality  (2.36),  the  equations 


y(yo)  =  y(T)  +  Z1 
Yl  =  T  +  Z2 


(3.83) 
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were  added  to  the  set  (3/63).,  Then,  Equation  (2.  35)  reads 


v  (T) 

11 11 

^12 

% 

V 

d,T 

0 

§;(T) 

^31 

^32 

”33 

”’34 

;  : 

dPlo 

d 

II 

>* 

C(T) 

”41 

^42 

”•43 

”44 

dp2 

0 

;  = 

0 

y(T) 

”21 

11 2  2 

r|23 

”24 

dp3 

o 

> 

i 

0 

0 

0 

0 

dp4 

o 

AZ2 

and  the  lbwer  two-by-two  matrix  of  the  product 

-I  o  *  -T 

(V»  \l/0)  V^J  (  7  ) 

xr2'  y  xT2 

is  the  desired  matrix  A4  .  In  the  other  method,  the  first  three  equations 
of  (3. 84)  were  solved  in  the  form 


dT 

_ 

v  (T) 

”u 

-1 

^12. 

-1 

”l3 

”l4 

dP3 

O 

dp1 

o 

= 

?(T) 

”31 

”32 

”33 

”34 

^4 

O 

dP  2 

o 

£ 

”41 

”42 

”43 

”44 

Let  D  be  the  matrix  of  (3.  85),  and  A^,  A 2  and  Ag  be,  respectively,  the 
upper  left  hand  (3  x  3),  upper  right-hand  (3  x  2)  and  lower  right-hand  (2  x  2) 
submatrices  of  (3.  82).  Then  the  desired  matrix  is 

D’  A1 D  +  D'A2  +  (D«A2)«  +  Ag  . 
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T".  -sq.  vm 


(3.84) 


(3.  85) 


1 


(3.  86) 


(3.87) 


Computer  Results  (C-6) 

Figures  3r5  through  3-8  represent  one  trajectory  obtained  using  the  modified 
Newton-Raphson  method  and  explicitly  completed  partial  derivatives.  (Com¬ 
pare  these  with  the  extremal  of  Figures  4-37  through  4t39,  which  was  used 
as  the  initial  guess  at  the  optimum.  )  The  control  function  (Figure  3-6)  starts 
on  the  -  6-degree  bound  (maximum  lift  condition)  and  stays  there  for  about 

65  seconds.  This  is  reasonable,  since  the  radiative  heating  rate,  propor- 

12  5  3 / 2 

tional  to  v  *  ,  could  become  quite  large  if  the  density  term  p  '  were  hot 

kept  small  by  the  control  maneuver.  As  velocity  begins  to  decrease  and  the 
< path  becomes  shallower  (Figure  3-5),  dissipation  of  energy  becomes  impor¬ 
tant;  therefore,  the  control  moves  toward  zero  degrees,  the  maximum  drag 
condition.  After  the  peak  heating  rate  (Figure  3-7),  a  little  more  lift  is 
called  for,  for  ranging  purposes.  Then  the  control  goes  to  the  +16-degree 
bound,  following  the  rule  of  thumb  that  the  convective  heating  load  is  lighter , 
the  faster  the  re-entry  is  accomplished.  The  peak  acceleration  is  9.  5  g's  com 
pared  to  the  10.  3  g's  of  the  original  extremal,  and  comparison  of  the  heating- 
rate  curves  of  Figures  3-7  and  4-39  shows  that  the  major  difference  is  near 
the  peaks  of  the  radiative  curves  (peak  value  is  188  BTU  for  Figure  4-39  and 
182  BTU  for -Figure  3-7).  The  value  of  the  optimal  criterion  was  about 
27,  500  BTU  for  the  extremal,  which  was  reduced  to  27,  334  BTU  by  the  opti¬ 
mization  method.  This  points  out  the  flatness  of  the  J  (y)  surface  for  this 
problem. 

Figure  3-8  is  included  to  show  the  angle  cp  of  Subsection  C-3.  It  can  be 
interpreted  as  the  unconstrained  value  of  the  control  function,  and,  over 
the  unconstrained  subarc,  u  and  cp  are  identical.  At  the  terminal  time, 
cp  goes  to  tt  because  of  the  necessary  condition  pg  (T)  =  0. 

There  is  some  doubt  that  the  trajectory  of  Figures  3-5  through  3-8  is  a 
relative  minimum,  although,  most  likely  it  is.  Both  test  matrices  of 
Subsection  C-5  are  indefinite,  but  numerical  problems  could  be  the  cause 
of  this.  The  determinant  of  the  matrix  in  (3,  84)  is  very  small  (on  the  order 
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of  10  ),  so  its  inverse,  used  in  (3,85),  cannot  be  known  very  accurately. 

In  the  matrix  (3.  87)  the  middle  two  terms  buck  the  outer  terms,  causing 
loss  of  from  two  to  four  significant  digits.  This  is  very  severe,  considering 
the  accuracy  of  the  Original  elements  and' the  number  of  matrix  multiplica¬ 
tions  required  to  arrive  at  the  result.  The  primary  cause  of  the  trouble 
appears  to  be  the  flatness  of  the  j  (y)  surface  for  this  problem.. 

Some  workers  in  the  field  feel;  that  the  initial  value  method  of  solving  the 
optimization  problem- is  too  sensitive  to  use  (for  an  example^  see  Ref.  29). 
They  argue  that  the  Euler -Lagrange  equations  are  the  adjoint  system  to  the 
original  system,  and  that  one  of  these  gives  unstable  solutions.  This  may 
well  be  true  for  some  problems,  but  no  such  instabilities  were  found  for 
this  application.  The  method  worked  well,  and  was  strongly  convergent  to- 
the  solution.  Adequate  prediction  for  all  variables  was  obtained;  in  fact,, 
the  desired  terminal  altitude  and  range  were  achieved  to  at  least  eight 

significant  digits.  The  terminal  value  of  p9  fared  somewhat  worse, 

-3  "  6 

7  x  10  compared  with  its  original  value  of  10  .  The  Hamiltonian,  a 

-7 

crude  measure  of  the  integration  accuracy,  changed  from  -3  x  10  to 
-2 

10  at  the>end  of  the  path,  an  acceptable  value,  according  to  past  experience. 
Other  evidence  supporting  the  numerical  accuracy  may  be  found  in  Table  3-1. 


Table  3-1.  A  Comparison  of  Partial  Derivatives 
Computed  in  Various  Ways 


- 

aj 

3J 

dJ 

3J 

ap1 

0 

o'P2 

0 

^3 

0 

&P4 

o 

Estimated 

-0.  14988xl05 

-0.  2i 7 x  10' 1 

0.  184x  10'3 

0.  118xl07 

Computed 

directly 

-0.  15057908  x  105 

0.  14247518  x  10_1 

0,  17744375  x  10"2 

0.  12049032x  107 

-p»(T)v  x(T) 
po 

-0.  150811xl05 

+  0.  142588  x  lO"1 

*  ' 

_ _ 

+  0.  177867x  10"2 

+  0.  12067  x  107 
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These  results  were  obtained  for  the  extremal  of  Figures  4-37  through  4-39.  The 
estimated  values  were  computed  for  program  checkout  purposes,  and  the  in¬ 
accuracy  of  the  second  result  comes  from  differencing  solutions  which  are  the 
same  to  seven  significant  digits  out  of  a  possible  eight.  Solutions  were  the  same 
to  six  significant  digits  for  ti  e  last  two  results,  but  to  only  three  for  the  first. 

The  last  row  was  handr computed  from  computer  results  truncated  to  eight 
significant  digits. 

Table  3-2  is  included  to  illustrate  the  convergence  of  the  modified  Newton^ 
Raphson  method  for  this  problem.  The  criterion  J(T)  rose  in  this  final  series  of 
iterations  because  the  trajectory  range  had  to  be  lengthened,  and  ranging  always 
increases  total  heat. 


Table  3-2.  Convergence  to  the  Solution 


Iteration 

number 

Method 

Constant 

Hamiltonian 
H(o)  . 

Criterion 

J(T) 

1  ■ . 

m$m Wm 

h(T)-A 

feet 

e-R 
;  feet 

1 

13.8 

26,  774 

1.  8  x  104 

0.  77 

-234,  691  . 

:  2 

0.  2 

11.  0 

I  26,887 

V,  5xl04 

0.  559 

-187,  674  ■ 

3 

0,28 

7.9 

27,019 

4 

1x10 

1.  121 

-132,  679  ' 

4 

0.46 

4.3 

27,  170 

5.  4x  103 

2.313 

-69,  331 

5 

0.  96 

0.  15 

27,  334 

-28.  0 

3.077 

-632 

6 

1.0 

-0,.2x  10"5 

27,  333 

27.0 

-0.  228 

-239 

< 

7 

1.  o 

0. 84  x 10"8 

27,  334 

0.  88x  10"1 

zero 

<  -0.  1 

8 

1.0 

-0.  65x  10"8 

27,  334 

0.  94 x  10"2  • 

zero 

zero 

9 

1.0 

-0.  27x  10'8 

27,  334 

0.  74xl0-2 

zero 

zero 
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SECTION  IV 

LINEAR  CONTROL  SYNTHESIS  AND  SIMULATION 


DISCUSSION  (A) 


It  is  well  known  that  a  quadratic  integral  criterion  when  applied  to  linear 
control  systems  specifies  a  linear  control  law  with,  in  general,  time 
dependent  gains.  Because  of  the  simplicity  Of  analyzing  and  mechanizing 
such  systems,  linear  control  in  the  vicinity  of  a  predetermined  reference 
trajectory  is  studied. 

If  the  equations  of  motion  for  state  x  under  control  u  are  represented  by 


I 


x  =  f  (x,  u)  (4.  1) 

and  a  control  u(t)  which  transfers  the  state  from  given  initial  conditions 
x  at  t  =  a  to  a  final  state  x. -  at  a  time  T  has  been  found,  the 
equation  for  the  correction  6  x  from  the  reference  path  x  (t)  to  a 
disturbed  trajectory  is 

6x  =  f(x  +  6x,  u  +  Au)  -  f(x,  u),  6x(a)  =  Axq,  (4.2) 

where  a  new  control  function  is  represented  as  u  +  Au.  Equation  (4.  2) 
is  the  exact  variational  equation.  In  the  neighborhood  of  the  reference, 
the  linear  approximation  to  Equation  (4.  2)  will  predict  a  correction  Ax 
which  should  reasonably  approximate  6  x.  The  "  first  variational  equation" 
is 


Ax  =  V  f  •  Ax  +  V  f  .  Au,  Ax  (a)  =  Ax  ,  (4.3) 

x  u  '  o 

where  the  matrices  of  partial  derivatives  V  f  and  V  f  are  evaluated 

r  xu 

along  the  reference.  The  problem  is  to  find  a  correction  to  the  control  Au 
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so  certain  final  conditions  for  the  reference  are  also  satisfied  by  the  perturbed 
trajectory.  The  simplest,  though  perhaps  not  the  most  natural,  condition  is 
to  require  particular  components  of  the  state  to  agree  at  the  end-point. 

Thus,  say 

A  x.  (T)  =  0  ;  i  =  q  +  1,  . , .  n  ;  0  =  q  =  n  r  1  .  (4.  4) 


Then  Au  will  be  completely  specified  as  a  function  of  time  if  it  is  choSen 
so 

fT 

J  (Au)2  dt  (4.5) 

a 

is  a  minimum,  Conditions  (4.4)  and  criterion  (4.  5)  can  be  greatly 
generalized,  and  the  derivation  of  the  control  law  may  be  made  using 
many  theories:  Pontryagin's  maximum  principle,  calculus  of  variations, 
dynamic  programming,  etc.  The  simplest  approach  for  the  case  at  hand 
will  be  used  here,  since  general  discussions  can  be  found  in  the  literature 
(Ref,  27  ,  for  example). 


The  region  of  validity  for  the  linear  approximation  is  found  by  experimentation. 
It  may  be  necessary  that  several  reference  trajectories  be  used  to  provide 
linear  control  over  the  entire  corridor.  This  depends  on  further  consider¬ 
ation  of  navigational  accuracy  and  model  fidelity.  In  a  subsequent  subsection, 
some  results  affecting  these  problems  are  discussed,  and  some  preliminary 
navigation  and  control  systems  are  evaluated  by  simulation. 


If  the  control  law  found  by  the  optimization  process  is  represented  by 


Au  =  E(t)  Ax  , 


(4.6) 


Equation  3)  for  the  linear  prediction  becomes  the  homogeneous  equation 


Ai  =  r 


V  f  +  V  f  .  E  Ax 
x  u 


(4.7) 
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If  all  state  components  must  agree  at  the  end-point  (q  =  0)>  that  point -must 
be  a  singular  point  of  the  differential  equation,  since  there  are  other  solutions 
than  the  trivial  Ax  (t)  =  6.  As  would  be  expected,  this  ease  ls  difficult  to 
control.  For  the  case  studied,  range  and  altitude  of  the  trajectories  are  rer 
quired  to  agree,  and  velocity  and  flight -path  angle  are  free.  This  demands 
less  of  the  control,  but  the  end-point  is  still  singular,  since  the  solutions 
begin  as  a  four-parameter  family  and  end  as  a  two-parameter  family.  This 
implies  that  some  of  the  gains,  the  components  of  E,  must  tend  to  infinity 
at  the  end-point,  which  will  lead  to  difficulties  in  simulation  and  mechani¬ 
zations.  Various  expedients  were  tried,  and  these  are  discussed  in  sub¬ 
sequent  subsections.  A  change  of  the  mode  of  control  near  the  end  is  probably 
indicated,  since,  in  any  case,  the  region  of  validity  for  the  linear  approxi¬ 
mation  is  exceeded  near  the  final  point. 


DERIVATION  OF  THE  CONTROL  LAW  (B) 


Let  the  row  vector  X  J  be  the  solution  of  the  adjoint  homogeneous  equation 

M  =  -X!  .  Vvf  ,  (4.8) 

J»  1  A 

which  has  all  components  zero  at  t  =  T  except  for  the  ith  component, 
which  is  unity.  By  direct  calculation,  using  (4.  8)  and  (4.  3),  it  is  found 
that 

d 

(M  •  Ax)  =  X!  V^f  .  Au  ;  (4.  9) 


and,  integrating  and  applying  the  end  conditions  (4.4),  the  equations 


X  J  (a)  .  Ax  (a) 


V  f  •  Au  dt 


0 


(4.  10) 


for  i 


q  +  1,  . . .  , 


n  result. 


The  control  must  minimize  integral  (4,  5),  subject  to  these  conditions..  This 
simple  problem  can  be  considered  using  constant  Lagrange  multipliers  p. 
and  completing  the  square  on  the  auxiliary  functional 


T  n 

J  (Au)2  dt  +  ^  p.  X!  (tQ)  Ax  (tQ) 
a  i=  q  +  1 


n 


T 


t  .  Pi 


i=  q+1 


Xf  •  V  f  .  Au  dt 
1  u 


a 


(4.  11) 


to  get 


-T  r 


a  L 


Au  + 


n 

T  £  pi  M  ’  V 

i  =  q  +  1 


dt 


+  l  terms  not  containing  Au 


(4.  12) 


The  terms  which  do  not  contain  Au  cannot  be  influenced.  Hence,  the 
minimum  occurs  for 

n 

4»!-t  [  pi H  •  V  •  <*• «> 

i=  q  +  1 
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Then  Equation  (4.  13)  reads 


Au  =  — p*  A‘  Vf  ,  (4, 16) 

and  (4.  14)  and  (4. 15)  become 

Y(T,  a)  p  =  2  A*  (a)  A  x  (a)  (4, 17) 

Y  (T,  a)  =  j  A'  y  (A’  •  y)'  dt  .  (4. 18) 

a 

The  matrix  Y  (T,  a)  is  assumed  to  be  nonsingular  for  a  <  T.  This 
property  is  related  to  the  controllability  of  the  system.  Thus,  the 
control  law  becomes 


Au-(t)  =  -  [  Vuf  (t)  ]  •  A  (t)  Y1  (T,  a)  A»  (a)  Ax  (a),  (4.19) 

when  the  transpose  of  Equation  (4.  16)  is  used.  This  holds  for  a  =  t  =  T 
and  may  be  considered  as  an  open-loop  control,  identifying  a  with  t 
gives  the  closed-loop  operation,  and  the  control  law  is 

Au  (t)  =  E  (t)  Ax  (t)  ,  (4,  20) 


with 


E(t)  =  -  [  vuf  (t))»  A  (t)  Y“!  (T,  t)  A«(ty  . 


(4.21) 
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SIMULATION  OF  THE  LOW-LIFT  VEHICLE  WITH  ROLL  MODULATION 
AND  ad  hoc  REFERENCE  TRAJECTORY  (C) 

The  Model  (C-l) 

The  motion  of  an  Apollo-type  capsule  is  assumed  to  be  governed  by  the 
two-dimensional  Equations  a^d  control  for  flights  at  a  fixed 
angle  of  attack  is  to  be  effected  by  changing  the  roll  angle  cp  .  It  is 
assumed  that  the  essential  characteristics  of  the  motion  are  preserved 
when  the  lateral  motion  and  forces  of  the  more  realistic  three-dimensional 
flight  are  neglected.  Hence,  the  lift  and  drag  coefficients  are  taken  in  the 
form  (1.  5).  The  atmosphere  is  generally  considered  to  be  exponential 
with  equation 

p  =  pQ  C  '  ,  (4. 22?) 

but  the  1959  ARDC  atmosphere  is  also  used  in  places  to  assess  variation 
in  density  effects. 

The  constants  in  the  equations  are: 

Aerodynamic  reference  area  : 

Weight : 

Drag  Coefficient:  ' 

(The  drag  coefficient  is  constant  for  Mach  numbers  greater  than  four;  it 
was  found  that  allowing  it  to  change  with  Mach  number  had  little  effect  on 
the  trajectory.  Hence,  it  is  taken  as  constant  throughout  the  analysis. ) 

Lift  coefficient:  --  CL0  =  0.  51  ^lq/Cqq  =  0,  5) 

Ballistic  coefficient :  W/C^  S  =  61.  7  lb/sq.  ft 


Sea  level  density: 
Atmosphere  coefficient: 
Sea  level  gravity: 

Radius  of  earth: 


p,  =  0..0023769 

p- 

ft  =  1/123500 

gQ  -  32.  2  ft /sec  ^ 

R  =  20;  903,  520  ft. 


The  initial  conditions  for  the  trajectories  are  chosen  as: 


Velocity:  vq 

Flight  path  angle :  7 

Altitude :  h 

o 

Range:  CQ 


The  terminal  conditions  expected  are: 

Altitude :  hT 

'-T  -  Range:  C 


36,  080  ft/ sec 

r  6.  4  degrees  (for  standard 
reference  trajectory) 

4.00,  000  ft. 

0 

100,  000  ft. 

2100  statute  miles 


The  time  of  flight  is  chosen  experimentally.  This  is  discussed  in  the 
next  paragraph , 

The  rate  of  roll  is  limited  to  20  degrees/ second  in  the  reference 
trajectories.  While  no  limit  is  placed  on  the  closed-loop  rate  +  A<p  , 
practical  considerations  require  that  this  limit  not  be  greatly  exceeded. 
This  condition  is  satisfied  on  the  flight  paths  studied.  The  control  com¬ 
mand  variable  A cp  is  restricted  by  the  condition  that  cp  +  Acp  must  fall 
in  the  same  half-circle,  either  ( o,  n)  or  (tt,  2 tt ) ,  as  cp  in  order  that 
an  artificial  control  reversal  does  not  occur  because  of  the  nonlinear 
manner  in  which  the  control  enters  the  equations  through  the  cosine 
function.  An  indirect  restriction  is  that  the  pilot’s  acceleration,  given  by 
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riot  exceed  10  g. 


The  ad  hoc  Reference  Trajectory  (C-2) 


To  begins  the  study  and  map  out  the  re-entry  corridor;  trajectories  were 
found  experimentally  by  guessing  forms  of  control  programs  q.'<t) . 
Trajectories  were  then  computed,  beginning  at  a  . height  of  400,  000  feet 
with  a  36,  080-feet  per  second  initial  velocity  and  various  initial  flight-path 
angles  yQ  .  It  was  determined  that  the  safe  re-entry  corridor  was  limited 
to  trajectories  beginning  with  y  between  -5.  4  and  -7.  4  degrees,  since, 
for  V  greater  than  -5,  4  degrees,  there  was  not  sufficient  lift  to  avoid  a 
skip-out,  and,  if  y  was  less  than  -7.4  degrees,  the  10-g  survival  de- 
acceleration  limit  was  exceeded. 

The  form  of  the  roll  program  adopted  is  as  follows:  The  roll  angle  cp  is 
set  at  zero  for  maximum  positive  lift  until  the  tangential  acceleration  at 
reaches  0.  6  g;  for  purposes  of  the  reference  program,  time  is  re-set  to 
zero  at  this  instant.  The  roll  is  held  at  zero,  for  15  seconds  past  the  re¬ 
set  time  and  then  is  increased  linearly  with  time  for  25  seconds  until  a 

certain  maximum,  cp  ,,  is  achieved.  This  maximum  is  held  for  the 

max 

next  80  seconds,  after  which  the  roll  is  decreased  linearly  for  80  seconds 

until  an  angle  of  22  degrees  is  reached.  It  is  kept  at  this  level  for  the 

remainder  of  the  flight,  which  is  terminated  when  the  altitude  decreases 

to  100,  000  feet.  The  value  of  cp  is  chosen  so  that  the  flight  has  a 

range  of  2100  miles.  A  plot  of  y  against  the  value  of  cp  ,  which  gives 

o  m  ctx 

the  2100-mile  range,  is  shown  in  Figure  4-1.  It  is  observed  from  this 


graph  that  a  constant-range  controller,  using,  the  -6,  4  degree  case  as  a  base, 
can  control  most  of  the  re-entry  corridor  (-7,  4  to  -5.  9)  with  a  maximum  Aty 
correction  of  less  than  15  degrees. 

The  roll  program,  for  the  trajectory  beginning  with  y  q  =  -6.  4  degrees 

is  shown  in  Figure  4-2;  this  is  taken  as  the  standard:  control-  program.  Its 
analytical  representation  is; 


Cp  (t) 

=  0,  0 

$  t  s 

15  ; 

cp<t) 

=  cp 

max 

t  -  1,5 
25 

,  15 

=  t  =  40  ; 

cp  (t) 

=  cp 
Tmax 

,  40  = 

t* 

120; 

(4. 

Cp  (t) 

=  cp 

max 

+  (22  - 

^max^ 

t'  120  120  1 
80  ’ 

t  =  200  ; 

cp  (t)  =  22  ,  200  =  t. 


The  quantity  cPmax  for  this  trajectory  is  90,  97  degrees,.  The  pilot’s 
acceleration  for  three  trajectories  with  rpll  programs  of  this  type  is 
shown  in  Figure  4-3,  These  trajectories  cover  the  entire  re-entry 
corridor  for  this  family  of  roll  programs.  The  re-set  time  is  plotted 
against  initial  angle  in  Figure  4-4. 

Characteristics  of  the  standard  reference  trajectory  (yq  -  6.  4  degrees) 
are  given  in  Figures  4-5,  4-6,  4-7,  4-8.  These  are  plots  at  height,  velocity 
and  re-entry  angle  against  range,  and  range  against  time  after  re-set.  The 
closed-loop  flights  described  later  use  this  as  a  reference. 


Sensitivity  of  the  standard  trajectory  to  various  types  of  errors  is  illustrated 
in  Figures  4-9  through  4-16.  Range  effects  from  shifting  the  roll  timing  ±2 
seconds  indicate  that  roll  timing  is  critical  (Figure  4-9).  Roll  magnitude 
variations  of  ±2  degrees  are  shown  in  Figure  4-10,  and  the  effect  of  a  change, 
of  ±0.  5  degrees  of  initial  re-entry  angle  on  range  is  given  m  Figure  4-Tf. 

Figure  4- 12  shows  the  effect  on  range  when  the  initial  height  on  the  standard 
reference  is  changed  by  either  +  20,  000  feet  or  -  50,  000  feet.  Skip-out 
occurs  with  an  altitude  error  much  above  20,  0,00  feet  in  an  open-loop  run. 

Initial  velocity  changes  of  ±1000  feet  per  second,  ballistic  coefficient 
variations  of  ±10  . per  cent  and  change  in  density  by  a  factor  of 

h  -  110i 000 

1.  0  ±0.5  - - - - —  ,  (4.  25) 

290,  000 

are  plotted  in  Figures  4-13  through  4-15.  Finally,  Figure  4-.16  shows 
the  range  effects  of  varying  the  time  that  cPmax  is  held  in  the  standard 
roll  program  (Figure  4-2)  by  ±20  seconds.  By  comparing  this  result  with 
the  runs  in  which  the  complete  roll  timing  program  was  changed  by  ±2  second's 
(Figure  4-9),  it  can  be  seen  that  the  range  effects  are  about  three  times  as 
great  as  those  in  Figure  4-16,  even  though  the  timing  change  in  Figure  4-1(5 
is  10  times  larger.  Roll  changes  applied  at  a  later  time  have  less  effect  on 
ranging  than  the  same  roll  changes  applied  earlier. 

Some  of  the  considerations  leading  to  this  choice  of  roll  programs  are  given 
in  the  following  discussion  : 

The  roll  cp  must  be  set  equal  to  zero  to  obtain  maximum  positive  lift 
during  the  initial  re-entry  phase,  at  least  until  the  first  acceleration  peak 
is  passed,  to  prevent  the  vehicle  diving  into  the  atmosphere  and  causing  an 
increase  in  magnitude  of  the  peak.  The  first  acceleration  peak  for  the  -  7.  4- 
degree  trajectory  is  equal  to  the  -  10-g  acceleration  limit. 
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After  the  first  acceleration  peak  is  passed,  the  vehicle  would  soon  skip  out 
of  the  atmosphere  if  the  roll  angle  were  left  zero.  Within  the  bound  of  the 
Apollo  roil-rate  limit  of  20  degrees  per  second,  it  was  necessary  to  apply 
rapidly  a  roll  sufficient  to  prevent  Skip-out. 

When  the  vehicle  slows  enough  that  a  skip- out  can  no  longer  occur,  the  roll 
must  be  decreased  to  a  small  value*  such  as  45  or  22  degrees,  to  prevent 
the  trajectory  from  going  too  deeply  into  the  atmosphere  and  exceeding  the 
acceleration  limit.  The  high  negative  lift  needed  to  prevent  a  skip-out  in 
the  case  in  which  y  equal  -5.4  degrees  cannot  be  left  in  effect  too  long,  or 
too  much  acceleration  will  occur.  The  final  value  of  roll  used  is  not  critical, 
because  most  of  the  ranging  is  done  by  the  maximum  roll  value,  it  is  best 
not  to  use  a  final  roll  of  zero  degrees,  because  no  additional  lift  is  then 
available  for  control  about  the  reference. 

By  basing  the  roll  control  on  the  interval  after  at  =  0.  6  g,  the  large 
variability  between  initial -time  and  the  time  when  the  vehicle  is  experiencing 
significant  aerodynamic  effects  is  eliminated  from  roll  control  timing.  With 
the  standard  roll  control  timing,  which  can  be  used  over  the  whole  corridor 
to  get  a  constant  range,  the  reference  is  less  sensitive  to  perturbations  and, 
also,  requires  smaller  A  <+>  corrections. 

This- roll  control  was  satisfactory  for  the  whole  corridor  if  timing  errors  are 
held!  to  within  a  few  seconds.  A  better  timing  would  start  the  roll  22  seconds 
after  the  a^.  =  0.  4  point.  This  would  delay  the  roll  application  in  the 

-  5.  4-degree  case  and  prevent  the  slight  increase  in  the  acceleration  peak 
of  the  aj.  =  0.  6  g  time  reset;  The  alternate  timing  is  set  to  leave  the 

-  5.  4-degree  case  unchanged,  because  any  delay  in  this  timing  of  more  than 
a  few  seconds  would  cause  a  skip-out.  Roll  timing  based  on  initial  time  is 
unsatisfactory  because  different  initial  conditions  would  require  initiation 

of  the  roll  program  to  occur  "over  a  range  of  40  seconds. 
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The  roll  angle  magnitude  which  will  maintain  a  constant-altitude  flight  at 
the  turnover  point  ( y  =  0)  is  an  important  quantity.  Increasing  this  value 
more  than  seven  degrees  will  cause  the  vehicle  to  dive  into  the  atmosphere 
so  deeply  that  the  10-g  limit  will  be  exceeded;  decreasing,  this  value  more 
than  seven. degrees  will  result  in  a  skiprout.  Ranges  between  2000  and 
5000  miles  are  possible  .for  all  initial1  conditions  by  using  different  values 
of  roll  within  seven  degrees  of  turnover  roll.  Ranging  possibilities  rapidly 
decrease  after  the  turnover  point,  especially  in  the  -  7.  4-degree  case.  It 
is,  therefore,  necessary  that  the  correct  roll  magnitude  be  closely  approxi¬ 
mated  in  the  vicinity  of  the  turnover  point  in  order  to  get  the  desired  range. 


The  Linear  Perturbation  Controller  (C-3) 

_  s  ,  ...  ' _  '  _ i _ -  _ / 


The  closed-loop  system  is  diagrammed  in  Figure  4- 17.  For  now,  it  is 
assumed  that  the  navigation  system  produces  an  accurate  estimate  xm(t) 
of  the  state  x(t).  (The  effect  Of  errors  in  navigation  is  discussed  under 
Subheading  C-6.  )  The  state  estimate  is  differenced  with  the  state  of  the 
reference  trajectory,  and  this  result  is  used  to  provide  a  correction  A  cp 
to  the  reference  roll  command.  The  inner-loop  controller  causes  the 
vehicle  to  assume  this  roll  angle;  this  controller  is  assumed  perfect  and 
is  not  studied  here.  The  related  equations  are 


m,,  v 

x  (t)  = 

x(t) 

Acp(t)  = 

E(t)  jjx^t)  -  xr(t)J‘ 

cp(t) 

<PC(t)  =  cpr(t)  +  Acp(t) 

(4.  26) 

dx 

dt 

f  <x,  cp(t)  ) 

r  p 

where  x  (t)  and  cp  (t)  are  the  reference  trajectory  state  and  corresponding 
control,  as  developed  in  the  previous  section.  The  feedback  gains,  repre¬ 
sented  by  matrix  E(t),  are  computed  according  to  the  theory  outlined  under 
Subsection  B.- 
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(The  quantities  p,  v,  y ,  C  are,  of  course,  those  of  the  reference;  the 
superscript  r  has  been  dropped  for  convenience.  )  The  guidance  sensitivity 
matrix  E(t)  is  calculated  according  to  formulas  (4.  18)  and  (4.  21)  from 
solutions  of  the  equations  adjoint  to  (4.  27).  These  components  are  plotted 
in  Figures  4-18  and  4- 19i 

The  elements  of  E(t),  from  the  theory  and  as  demonstrated  in  the  graphs, 
become  infinite  at  the  end-point,  and  it  is  necessary  to  introduce  some 
modification  to  take  care  of  this  singularity.  A  transformation, 

A  cp  =  arc  ;tan  A  cp  ,  (4.  28) 

was  used  for  this  purpose.  Experiments  indicate  that  this  transformation 
actually  increases  the  region  of  controllability.  The  additional  restriction 
discussed  under  C-  1  now  reads:  cp  +  A  cpc  is  limited  to  (0,  n  )  or  ( n,  2n), 
The  full  A  cp  control  can  be  used  during  the  time  cp  is  being  applied 
because  cp  is  close  to  90  degrees,  and  the  ±90  degree  range  of  Acp  will 
not  cause  a  crossing  of  either  the  zero-or  180-degree  limits  on  the  roll  angle. 


Closed-Loop  Results,  Precise  Navigation  Assumed  (C-4) 


Satisfactory  performance  was  obtained  over  almost  the  entire  corridor. 
Typical  terminal  errors  are  summarized  in  Table  4-1. 


Table  4-1.  Terminal  Errors  Caused  by  Various  Errors 
in  the  Model 


- 

Resulting  Terminal  Error 

Model  Variation 

Av(T) 

ieet/second 

Ah(T) 

feet 

AC 

miles 

1. 

cp  bias  -  1  degree 

15 

-21 

0 

2, 

’ft  i  a  -c  >h  ?  110,  000  \ 

■p'*  pl1  +  0-5  Wooo  ) 

173 

1791 

‘ 

0.  1 

3. 

'p  '*  p  K  1  '  °‘  5  290,  000  j 

-169 

-2584 

% 

i 

o 

• 

4. 

7q  =  -  6.  9  degrees 

86 

-  462 

0 

5. 

y f  =  -  7.  4  degrees 

180 

1282 

-°.  1 

6. 

vq  +  1000  feet/second 

-23 

-  53 

o 

• 

o 

< 

7. 

vq  -  1000  feet/sep;ond 

31 

-  49 

0 

8. 

1.  i  w/cD  S* 

-286 

- 10,  5  10r 

i 

o 

• 

to 

•  9. 

0.  9W/CdS* 

137 

4845 

-0.  5 

*The  constraint  on  Acpc  would  have  reduced  these  errors 


The  particular  case  in  which  the  initial  re-entry  angle  is  -  5.  9  degrees 
(compared  to  -  6.  4  for  the  reference  path)  and  in  which  the  exponential 
atmosphere  model  is  replaced  by  the  ARDC  '59  form  was  chosen  to 
illustrate  the  results.  The  roll  correction  Acpc  for  this  flight  is  shown 
in  Figure  4-20;  note  the  maximum  positive  lift  condition  of  -  22  degrees 
correction  after  370  seconds.  The  final  state  errors  are  increased 
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several  times  if  this  limit  is  not  imposed.  Figures  4-21  through  4-25  plot  the 
deviations  from  the  reference  for  range,  flight-path  angle*,  altitude,  velocity 
and  tangential  acceleration.  This  case  will  subsequently  be  referred  to  as 
the  standard  perturbed  trajectory; 

A  flight  beginning  with  angle  -  6.  9  degrees  is  graphed  in  Figure  4-26  and 
4-27.  Note  that  the  roll  corrections, are  more  in  the  linear  range  than  in 
the  previous  case  because  this  perturbation  is  toward  the  acceleration  limit 
boundary,  and  the  other  is=  toward  the  skip -out  limit. 

Figure  4-28  is  a  plot  of  AQ  ,  the  range  deviation,  against  time  for  two 
closed-loop  trajectories  with  the  same  conditions  as  the  reference,  except 
hQ  is  changed  by  ±50,  000  feet.  Since,  in  the  open-loop  case  with  ,Aho  =  25,  000 
feet,  skip-out.  occurs;  this  run  shows  that  the  linear  control  scheme  can  control 
large  perturbations. 

Other  perturbations  were  investigated.  The  results  are  summarized  in 
Tables  4-2,  4-3 , 


Table  4-2.  Perturbation  Study  Results  (1) 


Case 

AV(T) 

Ah(T) 

AC 

■I 

126 

688 

-  1,5 

2 

-197 

-2648 

0.  1  i 

3 

-330 

-6295 

-  4.0 

4 

-529 

-17885  ; 

-  2.  0 
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Case  1  is  the  standard  perturbed  reference.  Cases  2,  3  and  4  have  the  same 
initial  conditions  as  Case  1,  y  =  r  5,  9  degrees,  except  that  they  all  use 
the  exponential  atmosphere  instead  of  the  ARDC  *59  atmosphere.  Gases  1 
and  2  indicate  that  the  effects  of  using  the  exponential  atmosphere  ar.e  not 
large.  In  Cases  1  and  .2,  a  roll  correction  greater  in  magnitude  than  -  22  degrees 
is  called  for  after  the  time  the  reference  roll',  cp  is  set  to  its  final  value  of 
22  degrees.  Acp  is  then  set  equal  to  -  22  degrees*  so  that  cp.(t)  +  Acp  (t)  =  0, 
the  maximum  lift  condition.  In  Case  4,  Acpc  is  allowed -to  cross  the  zero 
point,  and;  in  Case  3,  Acp  is  set  zero.  Because  Case  2  has  final-  errors 
several-times  smaller  than  Cases  3  and  4,  it  is  concluded  that  the  best 
choice  is  to  hold  Acpc  at  -  22  degrees. 


Table  4-3.  Perturbation  Study  Results  (2) 


Case 

- 

AV(T) 

Ah(T) 

ARange 

5 

-  891 

-20,  508- 

-  26 

J6 

255 

455  . 

0.2  . 

7 

'  -1634  ; 

-63, 055 

-65- 

Gases  5  and  6  have  the  standard  reference  conditions,  except  the  initial  height 
is.perturbed  50,  000  feet,  in  Case  5  and  -  50,  000  feet  in  Case  6.  Case  5  shov/s 
the  effect  of  a  large  perturbation  towards  the  skip-out  boundary,  and  Case  6 
a  larger  perturbation  toward  the  acceleration  boundary.  Without  the  trans¬ 
formation  (4.  28),  it  was  found,  in  other  runs,  that  too  much  acceleration  for 
survival  occurs  (Case  5)  or  a  skip-out  occurs  (Case  6).  Both  Case  5  and 
Case  7  have  the  same  initial  conditions,  initial  Ah  =  50,  000  feet,  but  Case  7 

uses  the  unlimited  Acp  ,  In  Case  4,  the  use  of  the  unlimited  Acp  results  in  a 
several-fold  increase  in  final  errors,  but,  in  Case  7,  the  same  several-fold 
increase  of  the  already  large  final  errors  of  Case  5  may  affect  survival. 


A,  Perturbational . Navigation  System  (C-5) 


The  navigation  calculation  may  he  incorporated  into  the  control  scheme, as 
diagrammed  in  Figure  4-29,.  The  estimate  of  the  perturbation  of  the  state 
is  calculated  by  subtracting  .the  reference -tfafectory  acceleration  from 
that  which  is  measured  and  then  carrying  out  navigation  .computations. 

The  set  of  equations  to  estimate  the  state;pe4?t^^^U^4§i'thty§am & a| 

system  (4,  27)  when  the  A<|»  -term  ih  the  second  equation  of  that  set  is  - 

omitted.  Since  the  aerodynamic  force  terms  are  deter mihed  witb 
accelerometers,  a  more  convenient  form  for  the  first  two  of  these  equations 


•  Sn  :■  2  sin  y 

Av  =  -  Aa,  +  — -  •••  -k  :■  -  Ay  cos  y  +  — --■■■•  A  ?  : 
<1  +  0  :  '  1  +  5 


*  Aa 
Ar  =  H 


g  2  cos  y 

+  - - -w  :  Ay  sin  y  +  — - —  A  5 

v‘<l  +  5)  1  +  5 


g  cos  y  v  sin  y 

+  — - - -  Av  -  - - - - 

v  -(1+0  R  (1  +  0 


(4,29) 


cos  y  Ay 


1  +  5  R 


R  <1  +  5 ) 


A  mechanization  could  be  achieved  using  these  equations.  A  diagram  of 
such  a  system  is  given  in  Figure  4-30.  However,  not  all  terms  are  of 
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equal  -importance.  It  was  found'  that  the  equations  could  be  simplified 
to  the  form* 

Ay  =  -  Aa^, 

Ay  =  -i  i  Aa  -  (a  +  Aa  )  — 
v  |_  n  n  n  v 

,  •  -  i 

a?  =  (v  +  Av)  (y  +  Ay)  ^ 

-  1  1 
A  C  -  T+T  Av  * 

Note' that  some  second-order  terms  are  retained;  while  certain  first-order 
terms  are  dropped.  A  mechanization  for  these  equations  is  .shown  in 
Figure  4-31. 

'Using  the  approximate  scheme,  an  exponential  atmosphere  and  the  standard 
perturbed  reference  trajectory,  the  following  terminal  errors  were  found : 

Range  =  6  miles 

Altitude  =  '13,  02 5‘ feet 

Velocity  =  525  feet,  per  second  , 

A  further  analysis  shows  that  must  be  known  to  about  0.  2  degrees.  Since 
cannot  be  inferred  adequately  from  inertial  measurements  alone,  external 
information  is  required. 

The  accuracy  required  in  measuring  acceleration  is  indicated  by  inserting 
an  0.  05-g  bias  into  Aa,  and  Aa  .  The  results  are  given  in  Table  4-4. 

w  11 
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T+T 


E* 


go 


V  (1+?) 
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Av 


<4 
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Table  4-4,  Terminal  Errors  Caused  by  Accelerometer  Bias 
for  Standard  Perturbed  Trajectory 


Error  Source 

Terminal  Errors 

’  -c  •/ 

Actual. 

System 

Range 

•■r  V  C'  ' 

Altitude 

Velocity 

Range 

Altitude , 

Velocity 

Standard  perturbed 
trajectory 

,6.0  1 
miles 

-13,025 

feet 

;  -525 
.  ft/  second 

4  , 

:  Aa  =  Aa.  =  0. 05  g 
n  t 

54.5 

1000  , 

44 

-1,0 

miles 

,  -6000 
feet 

-920 

ft/ second 

:  Aa  =  Aa,  =  -  0.  05.g 
n  t 

o 
*  - 

o 

-20, 000 

-650 

_ _ _ J 

8.6 

-35,. 000 

,  -205 

A  control  system  which  uses  only  measurements  of  normal  and  tangential  accelera¬ 
tion  components  plus  a  roll  angle  determination^  was  developed  and  simulated.  It 
was  found  that  it  could  control  only  about  one -fifth  of  the  re-entry  corridor  and, 
so,  was  unsatisfactory. 


Effect  of;  State-Measurement  Errors  tC-6) 


The.  initial  conditions  for  the  flight  path  may  be  obtained  on  board  the  vehicle  from 
several  sources: 

Propagation  of  known  mid-course  errors 

Information  telemetered  from  a  ground  fix  near  the  onset  of  re-entry 
Instruments  in  the  vehicle. 

The  accuracy  and  availability  of  the  first  two  sources  for  estimating  the  initial 
state  are  not  known,  and:  it  is  doubtful  that  the  last  one  can  provide  the  flight 
path  angle  to  the  accuracy  required  for  successful  re-entry,  since  simple 


flight  path  indicators  generally  have  errors  of  about  one  degree.  The  Simula 
tions  reported  here  indicate  that  this  information  must  be  known  to  ±0.  2 
degrees  or  the  linear  control  scheme  may  fail,  initial  conditions  on  the 
state  other  than  the  flight  path  angle  are  not  critical.  Velocity  will  be 
known  to  100  feet  per  second  from  energy  considerations.  The  range  of  a 
perturbed5  trajectory  at  re-set  may  be  estimated  closely  by  multiplying  the 
range  of  the  reference  path  by  the  ratio  of  perturbed  time  to  re-set  to  the 
reference  time  of  reset  (re-set  time  is  the  instant  the  tangential  accelera¬ 
tion  reaches  0.  6:  g)„  From  this,  experience,  the  following  procedure  was 
adopted : 

For  the  control  calculations,  assume  initial  range  to  be  the  reference  range 
and  velocity  and:  altitude  at  re-set  time  to  be  that  of  the  reference  trajectory 
at  the  corresponding  instant.  The  range  at  re-set  is  estimated  as  outlined 
above. 

This  procedure  cauSes  errors  in  altitude  of  about  5,  000  to  10,  000  feet  becaus 
the  perturbed  p  is  based  on  the  ARDC  *59  atmosphere,  and  the  reference  p 
comes  from  the  exponential  atmosphere.  Since  the  ARDC  *59  atmosphere 
gives  a  -35  per  cent  to  20  per  cent  density  variation^ from  the  exponential 
atmosphere,  a  realistic  test  was  made  on  the  effect  of  density  variations 
on  the  controller. 

It  was  not  possible  to  handle  the  flight  path  angle  in  this  manner. 

At  one  stage  of  the  investigation,  it  was  hoped  to  estimate  the  perturbed  y 
when  a^.  =  u.  6  g  by  using  the  time  it  takes  a^  to  go  from  0.  2  g  to  0.  6  g. 

This  time  is  called  the  rise  time.  The  rise  times  for  a  series  of  initial 
re-entry  angles  were  determined  with  computer  runs.  Then,  by  an  inverse 
interpolation,  y  could  be  estimated  to  within  0,  15  degree  when  using  the 
exponential  atmosphere,  even  with  perturbations  such  as: 

Changes  in  the  ballistic  coefficient  of  ±10  per  cent 

Accelerometer  bias  less  than  0.  02  g 

Density  changes  by  a  factor  of  1  ±  0.  5  -  290^000^ 
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If  the  exponential  atmosphere  is  replaced  by  the  ARDC  *59  atmosphere,  an 
error  of  -1,  3  5  degrees  occurs  in  the  y  estimation.  Also,  it  was  noticed 

that  the  ARDG  ’82  atmosphere  gives  a  -0.  7  degree  error.  Tjoe  y  estimation 

* 

errors  .between  the  different  atmospheres  are  caused  by  changes  of  the  deriva¬ 
tive  of  the  density  with  respect  to  height.  Because  the  change  of  density  with 
respect  to  height  is  almost  inversely  proportional  to  molecular  temperature, 
and  this  temperature  fluctuates  a  great  deal  daily,  the  rise  time  of  a^  is  not  a 
satisfactory  estimator  of  y  .  When  the  atmosphere  is  changed  by  a  factor  as 
above,  the  height  at  which  ~  0.  6  g  changes  by  ±7000  feet;  but  a  similar 
change  in  the  height  occurs  at  the  point  at  which  a^  =  0.  2  g,  so  the  difference 

in  altitudes  is  still  within  800  feet  of  the  26,. 000-foot  difference  of  the  standard 
reference.  Using  thd  ARDC  *59  atmosphere,,  the  height  difference  is  19,000 
feet,  which  is  less  than  with  the.  exponential  atmosphere.  The  shorter  rise 
time  is  interpreted  as  a  steeper  flight  path. 

Many  closed-loop  trajectories  were  calculated  to  determine  the  effect  of  errors 
in  the  state  on  the  A^(t)  calculation.  The  standard  perturbed  closed-loop  tra¬ 
jectory  conditions  were  used  in  these  runs.  Table  4-5  gives  some  of  this  data. 


Table  4-5.  Effect  of  Errors  in  the  State  on  the  A$(t)  Calculation- 


Actual  Errors 

Estimated  Errors 

Case 

Av(T) 

Ah(T) 

AC 

Av(T) 

Ah(T)i  „ 

i 

126 

678 

-1.  5 

' 

- 

17 

-525 

-13025 

5.  7 

-702 

-20922 

2.1 

18 

44 

'  -1077 

54.  6- 

-920  , 

-6235 

-1.4 

19 

-851 

-19672 

-4P.  3 

-206 

-35306 

8. 0 

20 

-289 

-7982 

54.8 

-1072 

-32421 

8.5 

21 

-420 

-8902 

-42.1 

132 

-2939 

-0.4 

22 

25 

-1037 

-9.4 

3 

6603 

-11.6 

23 

-702 

-16048 

14.  8 

-939 

-39975  , 

1C,  3 
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Case  1  is  the  standard  perturbed  case,  which  uses  the  exact  state  in  calculating 
it  is  given  so  comparisons  may  be  made  with  Cases  17,  through  23;, 
which  are  under  the  same  conditions  except  the- approximate  state  error  is  used. 
The  integration  of  Equation  (4-2)  on  board  the  vehicle  gives  an  inexact  stale 
error  rather  than  the  exact  state  error  because  both  the  initial  state  error  is 
unknown  and  the  output  of  imperfect  accelerometers  must  be  used  in  calculating 
the  derivatives. 

Cases  XI,  22,  and  23  show  the  effects  of  integrating  Equaiioir  (4-2)  when  perfect 
accelerometers  are  used  but  the  initial  conditions  are  the  inexact  values.  In 
Case  17,  the  value  of  y  (0)  is  assumed  to  have  no  error;  in  Cases  22  and  23, 
the  error  in  the  y  (0)  estimation  is, assumed  to  be  0.  1  and  -0.  1  degrees, 
respectively.  In  Cases  18  and  19,  the  a^.  and  an  accelerometers  have  a  steady- 
state  bias  of  0.  05  g  and  -0.  05  g,  respectively.  In  Cases  20  and  21,  the 
accelerometers  have  2  per  cent  and  -2  per  cent  error.,-  respectively.  The 
initial  state  errors  in  Cases  18  through  21  are  estimated  in  the  same  way  as 
in  Case  17^  The  reset  times  in  the  last  cases  are  slightly  different  than  in 
Case  17  because  an  imperfect  accelerometer  is  being  used  to  determine  the 
point  at  which  a^.  =  0..6  g., . 

Integration  with  an  accelerometer  bias  of  0.05  g  for  600  seconds  gives  rise 
to  about  a  900-foot-per- second  velocity  error  and  a  50-mile  range  error, 
which  explains  the  50-mile  final  difference  in  A£  and  Av  in  Gases  18  and  19. 

The  controller  based  on  the  estimated  initial  state  does  a  good  job  of  reducing 
the  final  value  of  estimated  A£  to  around  five  miles,  but  it  is  off  about  50 
miles  from  the  real  range  error. 

Table  4.  6.  Perturbation  Study  Results  (4) 


Case 

Av(T) 

Ah(T) 

_ 

24 

25 

-549 

-733 

-11890 

-15623 

66.  5 

-78.  1 
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Cases  24  and  25  correspond  to  the  standard  perturbed  Case  1,  except  that, 
when  computing  the  control  A  $  in  Case  24,  Ay  was  omitted,  arid-,  in  Case  25,, 

Ah  was  omitted.  In  the  above  cases,  the  exact  state  is  known  but  not  complelely 
used.  When  both  Ay  and  Ah  are  ignored,  control  is  lost,  and  the  trajectory 
hits  the  ground  at  440  seconds  with  a  350-mile  range  error;  showing  the  neces¬ 
sity  of  estimating  the  complete  state  error.  If  the  perturbation  errors  get  too 
large,  for  instance  integrating  from  a  poor  path  angle  estimate,  the  controller 
fails; 

Figures  4-32  to  4-36  present  the  different  state  errors  in  Cases  1  and  20. 

These  are  plots  of  A  A£,  A  a  ,.  Ah  and  Av  with  time.  The  quantities  marked 
with  the  tilda  are  the  errors  computed  with  the  inexact  initial  conditions. 

The  difference  between  the  exact  and  inexact  error;  shown  by  Figures  4-34  and 
4-35,  is  small  for  about  the  first  100  seconds,  so  the  trajectory  is  controlled 
correctly  during  this  critical  time.  After  this,  the  approximate  state  estimate 
deteriorates,  but  bad  control  A^  (t)  at  this  time,  has  little  effect,  and  a  good 
re-entry  is  still  accomplished.  When  the  initial  Ay  estimate  is  in  too-great 
error,  ±0.4  degree,  the  altitude  estimate  rapidly  deteriorates,  and  erroneous 
control  signals  are  given  during  the  early  critical  part  of  the  re-entry.  This; 
sometimes  causes  the  trajectory  either  to  exceed  the  acceleration  limit  of 
10  g's  or  to  skip  a  thousand1  miles  or  more. 


STUDY  OF  A  LOW-LIFT  VEHICLE  WITH  ANGLE  OF  ATTACK  MODULATION 
AND  EXTREMAL  REFERENCE  TRAJECTORY  (D) 


The  Model  and  the  Controller  (D-li) 


The  re-entry  body  considered  in  this  subsection  is  the  capsule  studied  in 
Section  III  (C).  The  equations  of  motion  are  assumed  to  be  the  system  (1. 1), 
and  control  is  produced  by  changing  angle  of  attack  to  vary  lift  and  drag  coef- 
ficents  according  to  Equations  {1.  4>.;  The  numerical  values  for  the  coefficients 


in  these  equations:  are 


Cm  =  °'  88,  CDL  =  °‘  52'  cLO  =  -0*  505,  e  =  2,4. 

The  atmosphere  is. again  assumed  to  be  represented  by  the  exponential  low. 

A  linear  controller  is  studied,  the  theory  of  which  is  found  in  subsection  B, 

It  is  the  same  type  of  linear,,  closed-loop  controller  used  in  the  previous 
simulation.  The  control  attempts  to  drive  aititude  and  range  errors  to  zero, 
but  velocity  and  flight  path  angle  are  left  free.  Again,  because  of  the  singularity 
at  the  end r point,  special  procedures  are  employed1  to  limit  the  control  command 
in  that  neighborhood.  Two  methods  are  compared:  One  limits  the  control 
correction  predicted  by  Au  =  E(t)Ax  by  ±15  degrees;  the  other  bounds  the 
Au  near  the.  end  of  the  trajectory  by  1  1/ 2-times  the  maximum  correction  used 
in  the  beginning  of  the  flight.  (The  maximum  usually  occurs  at  the  first  turn¬ 
over  after  about  90  seconds  of  flight).  The  second  technique,  which  allows 
smaller-  variations,  was  generally  found  to  give  better  end-point  results  when 
the  disturbances  were  small;  but  with  large  perturbations,  the  results  are 
inconclusive. 


The  Extremal  Reference  Trajectory  (D- 2) 


The  linear  controller  is  intended  to  operate  about  the  optimal  reference  path 
which  was  computed  in  Section  III.  That  calculation  had  not  been  completed 
when  this  controller  study  was  begun,  so  ah  extremal  trajectory  was  chosen 
as  the  reference  for  an  intermediate  analysis. 

An  extremal  was  defined  in  Section  II  as  a  path  which  satisfies  the  Euler- 
Lagrange  equations  and  the  Weierstrass-Glebsch  conditions,  but  not  neces¬ 
sarily  the  required  boundary  conditions.  Here,  this  may  be  interpreted  as  a 
path  satisfying  the  equations  of  motion  with  the  control  required  by  the 


94 


minimum  principle.  As  in  Section  III,  the  control  command  was  limited  to 
excursions  between  ±16  degrees  (in  future  work  this  bound  is  to  be  extended). 


The  control  for  the  reference  is  plotted  in  Figure  4-37  along  with  the  pilot’s 
acceleration.  Figure  4-38  shows  the  state  variables  v„  a  ,  h  as  functions  of 
time,  and  the  heating  characteristics  for  a  re-entry  along  this  path  are 
recorded  in  Figure  4-39;  Note  the  similarity  of  these  graphs  to  those  for  the 
optimal  trajectory  (Figures  3-5  to  3-7). 

The  reference  was  computed  by  beginning  with  the  conditions 

vq  =  35,  000  feet  per  second 

y  =  -5;  75  degrees 
o  * 

hu  =  400, 000  feet 

O' 

C  =  zero  , 

*=0  ? 

and  a  particular  choice  of  the  adjoint  variables,  in  this  case 

=  3.  6565 
=  230 
=  67931.9 
=  -2565.  54  , 

and  terminating  when  the  second  adjoint  variable  became  zero.  The  range  of 
initial  guesses  on  the  adjoints  which  gives  reasonableire-entry  paths  was  found 
by  trial  and  error. 

The  control  sensitivity  along  the  reference  path  was  studied  by  applying  a 
constant  correction  to  the  angle  of  attack  (Au  =  e  Aa  =  2.4  Aa  )  and  in¬ 

tegrating  the  first  variai  -.ial Equations  (4.  3)  backward  in  time  from  the  end 


10 


20 


30 


40 


of  the  trajectory.  Typical  trends  are  shown  in  Figures  4-40  and  4-41  for  a  one- 
degree  change  in  angle  of  attack.  These  curves  may  be  interpreted  as  the  errors 
in  Av,  Ah,  arid  A£,  which  will  be  nulled  by  a  one-degree  correction  in  control' 
applied  at  the  time  indicated  by  the  abscissa  and  held  to  the  end.  Note  that  not 
much  control  may  be  exercised  near  the  end,  and  the  controllability  diminishes 
rapidly  after  the  first  90  seconds  of  flight. 


Open-  and  Closed-Loop  Results'  (D-3) 


The  controller  was  tested  by  calculating  a  series  of  trajectories  in  which  a 
single  perturbation  was  inade  in-an  initial  condition,  or  in  one  of  the  parameters 
of  the  vehicle.  Each  quantity  was  varied  by  small  positive  and  negative  incre¬ 
ments  and  by  large  changes  which  were  chosen  to  .be  2  l/'2-times  the  small 
variations.  The  open-loop  paths,  which  use  the  reference  program  for  angle  of 
attack  with  no  added  coritrol  correctiori,  were  coinputed  for  comparison. 

A  typical  case  is  illustrated  by  Figur  es  4r42  to  4-44;  the  initial  flight  path  angle 
was  increased  by  0.  05  degree  in  this  example.  The  control  correction,  plotted 
as  a  function  of  time  in  Figure  4-44,  saturates  at  about  220.  seconds,  and  was 
limited  to  -15  degrees,,  as  shown  by  the  dashed  line,  or  taken  as  -3.  91-degrees, 
drawn  in  a  solid  curve.  This  last  bound  is  1  1/ 2-times  the  control  correction 
peak  value  which  occurred  at  87.  seconds  flight  time.  The  final  altitude  error, 
shown  in  Figure  4-42,  was  reduced  from  568  feet  to  25  feet,  and  the  final  range 
error;.  Figure  4-43,  was  changed  from  0.096  mile  to  0,006  mile,  wheri  the 
smaller  bound  was  used.  These  numbers  are  compared  to  an  altitude  error  of 
2784  feet  and  a  range  error  of  22.  7  miles  for  the  open-loop  trajectory,  also 
illustrated  in  the  figures. 

Terminal  errors  resulting  from  perturbations  of  the  initial  state  are  listed 
in  Table  4-7.  For  each  change  two  closed-loop  runs  with  the  different  final 
control  bounds  and  the  openrloop  run  are  given.  Note  that  the  controller  allows 
altitude  errors  less  than  100  feet  and  range  errors  less  than  0.  01  mile,  usually, 
with  small  perturbations.  The  large  changes  are  controlled  within  a  few  thousand 
feet  in  altitude  and  a  half-mile  in  range. 
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Table  4-7.  Initial-State  Error  Cases 


>  Type  of 

'  Perturbation 

Maximum  A  u 
;  Correction^  in 
Degrees  . 

Final  Perturbed -State  Differences  from  Reference 

AV 

> 

ft/ sec. 

Ay 

degrees 

Ah 

'  feet  • 

;  AC 

miles 

A 

*c 

BTU’s 

A 

qr 

BTU's  ; 

None 

0.  33 

0.015 

-0. 003 

0.018 

0.0000 

0.01 

0.002 

- 

15 

^3.  765 

0.299 

1.  595 

s 

;  0,  0004 

0.01 

0.  002 

zero 

-0.036 

0.  002 

-1.205 

-0. 0048 

-0.05 

0.  005 

A  V  =  200  ft/  sec 

3.23 

J 

-11 

0.008 

47 

j 

0.008 

197 

442 

• 

:  15 

-11 

2 

485- 

0.081 

197 

442 

zero 

90 

4 

2602 

19.  892' 

535  ; 

392 

A  V  =-200  ft/  sec 

3.  21: 

4 

8 

r4 

-0.  001 

-193 

-417 

15 

zero 

2 

60 

0.  012 

-193 

0.  147 

. 

zero 

--67 

-4 

-2254 

-18.455. 

-509 

-377 

AV  =  500  ft/ sec 

8.  12 

-56 

-0.2 

-1051 

I 

o 

% 

h-L 

to 

, 

499 

1156 

15 

-57 

2.2 

-  606 

-0.  033 

499 

1156 

< 

zero 

295 

12.  0 

7378 

53.011 

1394 

1007 

AV  =  -500  fl/sec 

8.  00 

8 

2 

-14 

-0.  002 

-474 

-997 

15 

6 

3  .. 

21 

0.  006 

-474 

-997 

zero 

-142 

-9 

-5124 

-43. 864 

-1229 

-916 

97 
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Table  4-7.  Initial- State  Error  Cases  (Continued)' 


Type  of 
Perturbation 


Ah  =  -  2000  ft 


Ah  =  5000  ft 


-  137;  0. 020  1-179  -  96 


3320  33. 234  411  ,-192 


Ah  =  -  5000  ft 


5.  67 


-  16 


1.  2  -  956  -0.  130  157  110 


19  0.  006  157  110 


zero 


-  81 


4.9  -2763  -31.  563  -  401  183 


Ay  =  -0.  05  deg 


14  -  0. 002  29  103 


60  .  0.013 


29  103 


-2337  -20.  952  -  344  156 


Type  of 
Perturbation 


A*  =  0.05  deg 


&r  =  Oi,  15  deg 


Table  4-7.  Initial-State  Error  Gases  (Continued) 


....  .  Final  Perturbed -State  Differences  from  Reference 

Maximum  Au  "  '  • 

Correction  in 

s  Degrees 

-ft/sec  degrees  feet  miles  i  B.TU's  I  BTtPs 


3.  9: 
1! 

zert 


10.  27 
15 
zero 


25  0.006 

568  0.  096  . 

< 

2784  22.  715 


Ay  -j-0,.15  deg 

12.63 

o 

o 

-  2 

15 

-102 

-  1 

' 

zero 

-163 

-11 

- 

A  £  =  10  miles 

2.40 

-  24 

-  0.4 

15 

-  28 

3.  9 

zero 

zero 

zero 

A£  =  -  10  miles 

-2.40 

-  6 

■ 

0.  7 

15 

1 

2.5 

zero 

zero 

zero 

-  99 

-  99 
-156 


-0.  168  .  72  319 

-0.  115  73  319 

75.112  1151  -485 


2  -3013  — Q..3T8  -121  -277 

1  -3001  -0.  304  -121  -277 

IT  j  -6043  -58. 716  -954  437 


0.4  -  316  ,  -0.  039  -193  35 

3.  9  483  0.  090  -193  35 

ero  zero  10  zero  zero 


0.  7  -  413  -0. 049  200  -  33 

2.  5  54  0.  011  200  -  33 

sro  zero  -10  zero  zero 


Table  4  - 7.  Initial-State  .Error  Cases  (Continued) 


Type  of 

Perturbation  , 

Maximum  A  u 
Correction  in 
Degrees 

Final  Perturbed- 

State  Differences  from.  Reference 

AV 

ft/  sec 

■ 

Ay 

, 

degrees 

Ah 

feet 

aC 

miles 

A 

^c 

BTU'S: 

A 

<*r 

BTU’s 

AC  =  25  miles 

5.  98 

-118 

-5.5 

-4167 

-0.448 

-468 

88 

1 

1:5 

-131 

-3.  5 

-4481 

-0.459 

-468 

88 

zero 

zero 

zero: 

zero 

zero 

zero 

zero 

' 

A  C  =  -  25  miles 

6.  16 

-  21 

1.4 

-1338 

0.476 

515 

-81 

15 

17 

1,0 

.  39 

0.007 

513 

-81 

zero 

zero 

zero 

-  - 

.  zero 

-25 

zero 

zero 

The  first  three  cases  of  Table  4-7  are  a  check  on  the  accuracy  of  numerical 
integration,  since  the  reference  trajectory  and  control  are  computed  by  a 
'backwards  integration  from  terminal  time,  which  differs  from  the  forward 
calculation  by  small  amounts.  This  gives  rise  to  a  control'  correction  in 
the  closed-loop  case. 

Errors  produced  by  vehicle  parameter  changes  are  summarized  in  Table  4-8, 
The  controller  does  not  cope  aS  well  with  these  perturbations,  as  with  initial 
condition  errors,  because  the  control  correction  gains  Ett)  were' based  on.  the 
vehicle  parameters  of  the  reference.  The  final-state  errors  are  usually  less 
than  5,' 000  ft.  in  attitude  and  1  1/2  -  miles  in  range. 
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Table  4-8.  Vehicle  Parameter  Change  Cases  (Continued) 

.  Maximum  Au  Final  Perturbed-State  Differences 

iiion  Correction  i!1  AV  Ay  Ah  F  a!  Ta  Fa” 
Degrees  -  ^  qc  qr 

ft/ sec.  degrees  feet  '  miles  BTU  BTU 


10  per  centC 


-10*  per  centC^Q 


4  per  cent  C^jq 


»  -4  per  cent  CDO 


10  per  centC^Q 


-10  per  centC^Q 


4.  21 

17 

-  1. 7 

60 

-  0.  002  - 

4 

30- 

15 

9 

0.  3 

-  155 

-  0.  024  - 

4 

30 

zero 

zero 

-  67 

10.260 

183 

-  41 

5.  35 

-  21 

1.  5  - 

r  253 

-  0.  020 

8, 

34 

15 

,  -  24 

4.6 

'  321 

0;  073 

8 

34 

zero 

i“  7 

0.  1 

175 

-  9.798  - 

179- 

41 

1.05 

49 

2.3 

2222 

=  0.  205  - 

256 

!**  . . 

-  246  ; 

15 

76 

-  1. 4 

1722 

0.  202  ~ 

256 

-  246 

zero 

.-  13 

CD 

* 

o 

1 

275 

-  6.  292  - 

342 

-  24? 

1.06 

r  46- 

-  2.  3  , 

-  2267 

-  0.  144 

266 

261 

15 

-  78 

i 

0.2 

-3160 

-  0.  265 

206 

261 

zero 

14  , 

0.6 

-  280 

6.  542 

358 

256 

2,48 

182  , 

4.  9' 

5361 

0.  246  - 

624 

-  590 

15 

141  - 

1.8 

4419 

0.  511  - 

624 

-  590 

zero 

-  31 

1.5 

678 

-15.  298  - 

828 

-  580 

2.  51 

-  107 

5.7 

r  5595 

-  0.210 

683 

685 

15 

-  152 

4.  7 

-  7511 

-  0.554 

684 

685 

zero 

38 

1.6 

-  709 

16.  867 

926 

668 

Table  4-8.  Vehicle  Parameter  Change  Cases  {Continued) 


Final  Perturbed -State  Differences 


A 

BTU 

A 

BTU 

-  148= 

-  144 

-  148 

-  144 

-  198; 

—  141 

151 

149 

1 5  L 

149 

203= 

146 
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Comparison  of  Studies  in  Subsections  C  and  D  (D-4) 

The  studies  of  a  control  for  modulating  angle  of  attack  in  this  subsection  differ 
in  several  aspects  from  those  for  roll  angle  control  reported  in  Subsection  C. 

The  object  of  that  work  was  to  find  a  control  for  the  whole  re-entry  corridor, 
assuming  only  a  crude  navigational  system.  Very  large  perturbations  could  be 
controlled- t50,  000  feet  initial  altitude  error,,  for  example ---with  a  navigation 
calculation  from  biased  acceleration  data..  The  reference  path  was  chosen  so 
the  system  would  be  very  controllable. 

The  extremal  path,  on  the  other  hand,  was  found  to  yield  a  much  less  controllable 
system,  and  the  perturbations  had  to  be  strongly  restricted. 

MISCELLANEOUS  IDEAS  FOR  CONTROL  (E) 

The  linear  control  scheme  presented  appears  to.  provide  a  reasonable  engineering 
solution  for  the  re-entry  problem  considered.  There  are,  however,  a  number  of 
questions  which  must  be  answered  to  round  out  this  research  effort.  For  example, 
what  is  the  justification  in  using  the  linear  Equations  (4;  3)  for  prediction  so 
feedback  gains  may  be  found  when,  in  fact,  this  equation  holds  for  the  closed- 
loop  system  only  in  a  very  narrow  neighborhood  of  the  reference  and  not  in  ihe 
whole  controllable  corridor?  Why  is  improvement  noted  when  the  control  is 
modified  by  transformations  like  (4.  28),  or  bounds  imposed  when  these  modifica¬ 
tions  were  not  considered  in  the  basic  derivation  of  the  control  law?  What 
relations  are  implied  when  the  reference  trajectory  is  found  by  using  one  criterion 
and  the  feedback  control  found  by  using  another?  A  more  difficult  area  of  study, 
considered  a  little  in  Subsection  C,  is  to  find  the  minimum  information  of  the 
state,  number  of  variables  and  accuracy,  which  will  allow  a  successful  design. 

Some  of  the  answers  lie  in  a  more  careful  study  of  the  effects  of  the  essential 
nonlinearities  of  the  system.  A  step  along  these  lines  is  to  investigate  a  controller 
based  on  a  prediction  from  error  equations  linearized  in  the  state  variables, 
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but  in  which  the  nonlinear  aspects  of  the  control  authority  is  retained.  The 
control  criterion  must  also  be  changed  to  retain  the  nonlinear  character  of  this 
detail  of  the  system.  Or,  the  nonlinearities  may.be  avoided  by  adding  more 
reference  trajectories  and  error  equations.  If  two  paths  were  employed  as 
reference*  an  auxiliary  calculation  to  choose  the  reference  can  be  based  on  the 
construction  of  a  plane  midway  between  the  reference  states  and  perpendicular  to 
ihe  plane  formed  by  them. 

Another  approach  is  to  make  the  system  explicit,  calculating  a  new  optimal 
trajectory  at  each  point  along  the  motion.  This  may  be  done  by  simplifying  the 
Newton-Raphson  method  of  Section  II,  by  dropping  more  unessential  terms  using 
schemes  based  on  the  second  variation,  or,  by  using  a  Newton-Kantorovich-type 
of  iteration  to  predict  an  optimal  path.  Or,,  the  equations  of  motion  may  be 
solved  approximately  by  Galerkin’S  technique,  and  a  suboptimal  path  found  by 
the  Rayleigh-Ritz  procedure.  This  would  require  relatively  elementary  calcu¬ 
lations. 


SECTION  V 

CONCLUSIONS  AND  RECOMMENDATIONS 


ACCOMPLISHMENTS 

A  computational  procedure  is  developed  which  yields  the  control  program  and  the 
corresponding  re-entry  flight  path  to  minimize  a  given  criterion.  In  the  particular 
case  studied,  this  criterion  is  the  total  beat  to  the  vehicle.  The  calculation  is 
automatic,  in  that  the  calculation  will  run  on  a  large  digital  computer,  iterating 
as  need  be,  until  an  optimal  trajectory  is  found,  without  intervention  by  the 
operator.  The  basis  of  the  calculation  is  the  theory  of  the  calculus  of  variations. 

It  is  demonstrated  that,  when  suitably  modified,  this  theory  provides  the  flexi¬ 
bility  and  a  generality  to  encompass  a  large  class  of  control  synthesis  problems. 

The  chief  numerical  method  employed  is  a  modified  Newton-Raphson  technique. 

The  method  is  superior  to  both  the' gradient  approach. and  standard  Newton- 
Raphson  methods,  which  were  tried,  in  terms  of  rapidity  of  convergence,  ease  of 
automation  and  auxiliary  computations.  The  scheme  is  insensitive  to  round-off 
errors,  and  .it  exhibits  no  tendencies  toward, instability  because  of  adjoint  solutions. 
It  was  found  to  be  better,  for  accuracy  and  computer  usage,  to  compute  the 
necessary  partial  derivatives  from  explicit  formulas  derived  from  the  second 
variations,  rather  than  calculate  them  from  finite  differences. 

The  success  in  computing  optimal  trajectories  is  due,  in  part,  to  the  simple 
device  of  limiting  the  control  authority  so  that  a  reasonable  guess,  close  to  a 
ballistic  trajectory  and  in  the  region  of  convergence  for  the  method,  may  be 
made  and  subsequently  relaxing  this  limit  as  the  re-entry  corridor  is  studied. 

It  is  demonstrated  that  a  linear  perturbational  control  could  be  designed  to 
provide  control  within  a  substantial  re-entry  corridor  about  a  carefully  chosen 
reference  path*  The  accuracy  requirements  for  system  mechanization  are  not 


severe.  A  system  with  accelerometers  and  ah  rttitude  reference  of  reasonable 
accuracy  would  assure  safe  re-entries.  It  will  control  vehicles  with  errors  in 
■parameters  and  changes  in  atmosphere  density  profile.  Modification  of  the  linear 
controller  gains  at  the  end-point  is  examined  with  care. 

Several  complicated  digital  computer  programs,  'which  use  specially  developed 
integration  and  interpolation  methods,  were  constructed.  These  will  be  in¬ 
valuable  for  future  investigations. 

This  study  was  conducted  with  a  two-dimensional  model  of  the  re-entry  maneuver. 

This'  research  has  laid  the  theoretical  foundations  requisite  for  engineering 
feasibility  studies  of  optimal  re-entry  control.  Further,  limited  experience  with 
these  techniques  suggests  that  an  optimal  rer.entry  trajectory  has  real  and  significant 
performance  advantages  in  terms  qf  the  chosen  optimization  criterion.  This 
research  has  also  made  it  clear  that  on-board  computation  of  true  optimal 
trajectories  is  not  likely  to  -be  feasible  in  the  foreseeable  future.  However,  such 
trajectories  can  be  pre-computed  and  stored  on  board  for  use  in  the  iinear  control 
scheme.  Alternatively,  there  is  a  good  ^probability  that  simplifying  approxima¬ 
tions  can  be  made  to  the  theory  which  will  make  feasible  the  on-board  computation 
of  sub-optimal  paths.  The  merits  of  these  simplified  schemes  can  be  evaluated 
by  comparison  with  optimal  trajectories  computed  by  the  methods  developed  during 
this  research  program. 


RECOMMENDATIONS 

An  original  objective  of  this  program  was,  and  still  is,  to  examine  the  mechaniza¬ 
tion  problems  associated  with  optimal  or  near-optimal  re-entry  control  systems. 
The  prerequisite  theoretical  work  is  now  sufficiently  complete  that  these  problems 
may  be  studied.  The  program  suggested  below  has,  as  its  principal  purposes, 
the  development  of  a  technologically  feasible  mechanization  and  the  evaluation  of 
the  performance  degradation  which  this  sub -optimal  system  exhibits  with  respect 
to  an  optimal  system. 


Task  1;  Optimization 


There  are  still  some  aspects  of  the  Newton-Raphspn  optimization  scheme  which 
should  be  pursued.  These  are  minor  in  nature,  but  should  produce  dividends 
in  faster  determination  of  optimal  trajectories.  They  are: 

Examine  methods  for  increasing  rate  of  convergence  to  the  optimal 
trajectory 

Examine  other  differential  equation  integration  algorithms  for  faster 
path  solution. 

Examine  other  optimization  criteria 


Task  2:  Linear  Control 


The  various  linear  control  schemes  now  known  should  be  tested  against  each 
other  to  determine  their  effect  upon  the  optimization  criterion.  The  results  of 
these  studies  will  determine  the  linear  control  law  to  be  mechanized.  The 
linear  control  schemes  include: 

Minimum  deviation  of  control  from  reference  control 

Minimum  deviation  of  state  from  reference  state 

Optimal  control  in  the  vicinity  of  the  optimal  trajectory 


Task  3:  Approximation  Techniques 


The  mos.  fruitful  areas  of  on-going  research  will  be  those  of  finding  suitable 
simplifying  approximations  which  will  permit  near-optimum  trajectories  to  be 


computed  on  board  the  vehicle  and  of  evaluating  the  performance  of  the 
resulting  sub-optimal  control  scheme. 

Practical  utilization  of  the  theory  which  has  been  developed  depends  on  the 
successful  extension  of  these  results  toward  simplifying  approximations  which 
retain  as  many  performance  advantages  of  an  optimal  system  as  possible.  Some 
promising  approaches  are: 

Rayleigh-Ritz  techniques  for  optimization  arid  a  Galerkin  approximation 
to  solve  the  equations  of  motion 

Optirrial  solutions  of  simplified  equations  of  motion 

Djmamic  programming  techniques  for  on-board  computation  of  a  crude 
approximation  to  the  optimal  path 

Task  4-:.  Simulation 

A  re-entry  control  simulation,  designed  to  evaluate  the  most  promising 
techniques  resulting  from  Tasks  2  and  3  with  respect  to  engineering  feasibility, 
should  be  accorriplished. 
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APPENDIX.A 

THE  EXTREMAL  DIFFERENTIAL  EQUATIONS  AND  PARAMETERS 


! 
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An  extremal  is  defined  in  Section  II  as  a  path  satisfying  the  relations  (repro¬ 
duced  here  for  convenience)' 


■p»*  *  "  v  f  +  p»7  f  +  G 

(A.  1) 

0  =  v  f  +  p»v  1  +  fjL.‘v  G 

U  0  i  u  u 

(A.  2) 

o  =  fi.G.,  j  =  1, .....  q 
v  J 

(A.  3) 

0  ^  |JL 

(A.  4) 

X  =  f,  x(o)  =  xq 

(A.  5) 

•  2  ~  n 
a  =  G  £  0 

(A.  6) 

H(t,  x,  p,  u )  £  H(t,  x,  p,  U) 

(A.  7) 

TT*  V.  2  (H  +11G)  TT  I>  0 
u 

(A.  S’) 

V  G  n  =  0 
u 

(A.  9) 

It  is  convenient  to  consider  each  subarc  of  the  extremal  separately.  For  this 
purpose,  it  will  be  assumed  that  all  the  constraints  are  greater  than  zero  along 
the  first  subarc.  Then,  according  to  Valentine  (Reference  23)  the  inequality 
constraints  may  be  neglected:  over  this  subarc. 
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The  corresponding  Euler -Lagrange  Equations  (A.  1)  and  (A.  2)  are 


S V  f  +  p{7  f  =  V  F 

1  xo  r  X  X 


1 


(A.  TO) 


0  =  7  f  +  p'V  f  = 

no  e  u 


V  F 
u  1* 


(A,  11) 


It  is  well  known  from  implicit  function  theory  that  Equations  (A.  11)  can  be  solved 
(at  least  in  principle)  in  the  form. 


u  =  u(Xj  p) 


(a.  m 


if  the  determinant 
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(A.  13) 


is  different  from  zero.  Equation  (A.  13)  is  the  determinant  of  the  Hilbert 
differentiability  condition  for  this  problem.  Substitution  of  (A.  12)  into  (A.  10) 
and  (A.  5)  thus  produces  a  2n  set  of  differential  equations 


£  =  f^x^p),  x(0o)  =  XQ  (A.  14) 

4*  =  vxFl(x'p>*  (A.  15) 


These  will  be  called  the  reduced  differential  equ£  ions  of  the  extremals  for  thw 
first  subarc.  It  is  well  known  from  the  theory  of  differential  equations  that 
(A.  14)  and  (A.  15)  have  a  unique  solution  for  a  given  set  of  initial  conditions. 


:r 


But  half  of  these,  x0*  are  given  for  the  problem.  The  solution  is  thus  a  func¬ 
tion  of  the  n  initial  conditions  p  and:  is  said  to  be  imbedded  in  an  n-.parameter 

O 

family  of  extremals.  Over  the  first  subarc,  then,  the  solution  has  the  form 


x  ~  x(t,xo,  pQ)  =’  x(t,:pQ) 


P  =  p(t,  xq,  po)  =  p(t,  po) 


u  £  u(t,xo,po)  =  u(t,p0). 


(A.  16) 
(A,  17.) 
(A.  18) 


It  is  further  known  from,  the  theory  of  differential  equations  that  Equations  (A,  16) 

and  (A.  17)  and,  consequently,  (A.  18)„  have  continuous  partial  derivatives  in  the 

variables  !,  x  and  p  of  at  least  second  order; 
o  o 

The  above  results  are,  of  course,  obtainable  from  the  general  form  of  the 
problem  of  Bolza.  It  is  known  that  the  arc  can  beimbeddedln  an  (n+m)  parameter 
family  of  arcs,  and  that  there  are  2n  +  2m  differential  equations  of  the  extremals. 
The  2m  differential  equations  become  the  m  algebraic  Equations  (A.  11)  arid  their 
time  derivatives  (which  introduce  na.new  information)  when  initial  conditions  are 
imposed.  Furthermore,  the  determinant  in  the  Hilbert  differentiability  condition 
reduces  to  the  form  (A.  13). 


Over  the  second  subarc,  it  is  assumed  that  one,  of  the  inequality  constraints, 
say  Gj,  is  equal  to  zero.  The  function  (2.  10)  is  written,  following  Valentine 
(Reference  23),  as 


CO 

II 

fQ  +  jAf-x)  + 

(A,  19) 

and  the  Euler- 

Lagrange  equations  as 

-p.  . 

'Vo  +  e’V^iW  V2 

(A.  20) 

0  = 

Vo  +  P'V  +  Wl  =  VuF2- 

(A.  21) 

115 

To  this  is  added  the  equality  constraint 


0  •=?  Gr  (A,  22) 

If  the  m+i  Equations  (A.  21)  and  (A.  22)  are  to  be  solved  for  the  m+1  variables  u 
and  [4.^,  the  determinant 


det 


V2F9 

U’  2 


>v;  G. 
u  1 


7’G. 

u,  1 


(A.  23) 


must  not  be  zero.  Equation  (A.  23)  is  the  determinant  of  the  Hilbert  differenti¬ 
ability  condition  for  this  arc.  It  is  then  found  that 


u.  =  u(x,  p)  (A.  24) 

V-l=  (x,p)  (A.  25) 

and,  consequently,  that 

x  =  f2  (x,  p)  (A.  26) 

-p«  =  VxF2(x,p).  (A.  27) 


are  the  reduced  differential  equations  of  the  extremals  for  this  arc.  These 
equations  again  have  a  unique  solution  for  a  given  set  of  initial  conditions. 
Furthermore,  the  solutions  possess  continuous  partial  derivatives  of  at  least 
second  order  with  respect  to  t  and  the  initial  conditions  and  p^.  It  remains 
to  be  shown  that  the  solutions  over  this  arc  are  continuous  functions  of  the 
initial  conditions  pQ  and  that  partial  derivatives  of  at  least  second  order  exist. 
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The  terminal  point  of  the  first  subarc  is  defined  by  the  equation. 

Gl|  x  (tl»  P0)  *  P  (tr  9  o'*]  =  °*  (A.  28); 

This  may  be  solved  for  tj.(pQ)  if  the  derivative 
61  =  7xGl  *  +  7pGlP 

is  non-zero.  Equation  (A.  28)  represents  the  equation  for  the  terminal  surface 
of  the  family  of  extremals  which  are  ^solutions  of  Equations  (A.  14)  and  (A.  15). 

This  terminal  surface  represents  the  initial  surface  for  the  family  of  extremals 
over  the  second  subarc  which  are  solutions  of  Equations  (A.  26)  and  (A.  27). 

For  a  given., set  pQl  the  terminal  values  x  (t^,  po)i  p  (t  ,  pQ)  are  the  initial 
values  for  the  differential  Equations  (A.  26)  and  (A.  27).  This  follows  from  the 
continuity  of  x(t)  arid  pit),  Then,  since  the  solutions 


x  (t)  =  x  (t,  x1#  Pj  ),  t  *  t1  (A.  29) 

p  (t)  =  p(t,  Xx,  P1),  t  2:  t1  (A.  30) 

are  continuous  functions  of  x^  and  p^,  and  since 

xx  =  x  (ty  pQ)  =  x^p^  (A.  31) 

pl  =  p  <tr  pQ)  =  p!<p0>  (A,  32) 

are  continuous  functions  of  p  alone,  it  follows  that 

ro 

x(t)  =  x(t,  pQ),  t  ^  t1  *  (A.  33) 

p(t)  =  p(t,  pQ),  t  s  .  (A.  34) 
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Furthermore,  since  continuous  partial  derivatives  of  at  least  second  order 
exist  for  Equations  (A.  29)  -  (A,  32)  in  the  indicated  variables,  it  follows 
that  (A.  33)  and  (A.  34),  possess  partial  derivatives  of;at  least  second  order 
in  t,  p.i  Since  there  are  no  corners,  IA1  is,  at  least,  absolutely  continuous. 

Further  subarcs  may  now  be  added,  as  long  as  the  number  of  equality  constraints 
does  hot  exceed  m,  the  dimension  of  the  control  vector.  Each  subarc  possesses 
control  functions  and  multipliers  of  the  form  (A.  24)  and  (A.,25)  and,  furthermore, 
differential  equations  similar  to  (A.  26)  and  (A.  27).  It  is  readily  verified  that 
the  solutions  are  continuous  functions  of  the  initial  conditions*  pQ  and  t,  and 
that  continuous  par tiarderivatives  of  at. 'least  second  order  exist. 


Points  where  equality  constraints  change  to  the  "  greater- than"  state  will  now 
be  examined.  For  this  purpose  it  will  be  assumed,  that  only  one  of  the  con¬ 
straints,  say  Gr. ,  is  zero  over  the  first  subarc,  and  that  it  is  greater  than  zero 
over  the  second.  Equations  (A.  24)  -  (A,  27)  hold  over  the  first  .subarc,  and 
(A.  12),,  (A.  14)  and  (A.  15)  over  the  second  subarc.  The  Equation  (A,  28)  for 
the  terminal  surface  is  replaced  by  the  equation 


[x  <tr  Po),  P<t1,PQ)]  =  0, 


(A.  3  5) 


since  h-j  is  a  continuous  function  and  must  go  to  zero  before  can  be  greater 
than  zero.  It  is  then  seen  that  the  arguments  follow  through  as  before,  provided 
the  determinants  and  are  different  from  zero  and  that  the  derivative 


=  V  x  +  V  up 
x  1  p  1“ 


is  non- zero  on  the  terminal  surface  of  the  first  subarc.  On  the  second  subarc, 
then,  the  solution  is  a  continuous  function  of  t  and  pQ,  and  partial  derivatives  in 
these  variables  of  at  least  second  order  exist.  The  result  is  easily  generalized 
to  several  equality  constraints  going  to  the  "  greater-than"  state  over  a  series 
of  subarcs. 
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To  sum  up,  a, path  is  composed  of  a  finite  number  of  subarcs,  each  of  which 
has  its  own  set  of  differential  equations  and  terminal  surfaces.  The  first  subarc 
has  a  specified  set  of  initial  conditions,  xq,  and  the  last  has  a  specified  terminal 
surface.  Equations  (2. 4).  The  equations  of  this  terminal  surface  have  continuous 
partial  derivatives  of  at  least  second  order  with  respect  to  the  initial  condi¬ 
tions  p^,  and  the  terminal  value  of  the  independent  variable  T. 

Now  introduce  the.  new  differential  equation 


xn+,l  =  V 


x  ,  1  (  0  ),  =  0, 

n+1  o 


(A.  36)= 


where  f  is  the  integrand  of  Equation  (2. 1).  The  solution  to  this  is  easily 
0 

seen  to  be  of  the  form 


xn+l  =  xn+lU*  P0 


■(t,p  )  =  f  f  (x,p)dT 


AA.sn 


over  the  path  andi  on  the  terminal  surface. 


Vi1*  ■  x„+i<T-p0>-  (A-38) 

Furthermore,  since  g  of  Equation  ( 2,  1)  is  a  functional  of  continuous  functions, 

J(T)  =>  J(T,po)  =  g[T,x(T,po)]  +  xn+l(T’ po- •  (A-38) 

The  terminal  surface  Equations  (2,  4)  may  be  expressed  in  the  form 

\p(T)  =  i^<T,po)  =  ^[T,x(T,po)j  =  0.  (A.  40) 

Thus,  the  problem  reduces  to  that  of  minimizing  (A.  39)  subject  to  the  constraint 
Equations  (A.  40)  in  the  indicated  variables. 
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APPENDIX  B 

THE  SUFFICIENCY  CONDITION  FOR 
THE  FIXED  END-POINT  PROBLEM 


When  the  terminal  surface  consists  of  less  than  (n+1)  equations  of  the  form  (2.  4) 
it  is  always  possible  to  use  the  sufficiency  conditions  of  Section.  II  (C),  As  shown 
in  Section  II  (F)j,  the  second  derivative  matrix  comes  from  evaluation  of  the 
second  variation,  using  neighboring  extremals  (the  accessory  minimum;  problem) 
which  satisfy  the  necessary  conditions  on  the  problem*.  When  the  end-point 
is  fixed,  i.  e. ,  r  =  n+1,  it  is  in  general  impossible  to  satisfy  both  the  initial 
and  terminal  conditions  with  this  type  of  neighboring  extremal  (assuming  the 
problem  is  normal  and'-  that  the  x  matrix  is  non- singular  at  the  end-point). 
However,  it  is  still  possible  to  construct  neighboring  extremals  which  have 
discontinuous  derivatives  at  One  point  and  which  do  satisfy  end  conditions.  The 
second  variation  is  evaluated  with  these  neighboring  extrelhals  to  establish 
the  sufficiency  condition  for  a  relative  minimum  for  this  problem;  Bliss' 
theorem  86. 1  establishes  the  basis  for  the  development.  ** 


The  complete  fundamental  solution  matrix  of  the  system  (2.  80)  and  (2.  81)  is 


V  X 

xo 

V  x"l 
P0 

II 

nil 

i - 

CM 

V  p 

L  0 

V  p 
p0  . 

n21 

n22 

m 

where  the  tt  notation  is  introduced  for  the  convenience  of  the  following  develop¬ 
ment.  The  initial  condition  for  (B.  1)  is  the  2n  x  2n  identity  matrix.  Then  any 
solution  6x(t)  and  6p(t)  which  starts  with  6x(0)  =  0  must  be  composed  of  a 


*See  Bliss,  reference  22,  pp.  226-234,  243-247,  and  253-257  for  the  discussion 
presented  here. 

**Ibid,  p.  246 
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linear  combination  of  it  ^  ^ 

aadn22: 

^  6x(t)  =  Tr12^t>a, 

6x(o)  =  0 

(B.  2) 

5p<t):  r  ^2(tja, 

6jp(o)  -  a. 

(B.3) 

‘  U(T)" 

r'  ■ 

'"„ITI  "12«T) 

■*  - 

c 

V(T) 

tt21(T>  n22(T) 

D 

where  a  is  an  arbitrary  (nxl)  dimensional  vector.  It  is  desired  to  construct 
similar  neighboring  extremals  which  have  the  property  that  6x(T)  =  0.  This  is 

done  by  requiring 


(B.4) 


where  C  and  D  are  (nxn)  matrices  to  be  determined,  Z  is  the  zero  matrix  and 
U  and  V  are  matrices  such  that 

6x(t)  =  U(t)b,  6x(T)  =  0  (B„  5) 

,6p(t)  =  V(t)b,  6p(T)  =  b,  (B.  6) 

and  b  is  an  arbitrary  (nxl)  dimensional  vector.  Now  matrix  {B.  1)  is  non-singular 
by  definition,  so  it  may  be  inverted.  Call  the  inverse  at  t  =  T  the  A  matrix. 

Then,  from  (B.  4) 


*1 

C 

A11  A 12 

r  ^ 

z 

A12 

D  * 

_A21  A22 

I 

A22 

<B,  7} 


from  which  it  follows  that 


U(t)  =  nn<t)  A12  +  n12(t)  A22 

(B.  8) 

VM  =  n21<t)A12  +  TT22(t)A22. 

(B.  9) 
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Other  relationships  required  ini  the  development  are 


T72i'  11 14 

V 

^21 

TT22*  nl2 

TT12' 

n22 

TT22I  nll 

^12* 

n2l 

(B,  10) 
(B.  11) 

(B.  12) 


which  hold  at  any  point  t.  These  relationships  may  be  verified  by  differentiating 
and  substituting  Equations  (2. .80)  and  (2.  81)  into  the  result.  It  is  readily  verified 
that  Oj  =  Cg  =  Z  ahd  that  =  I  for  the  initial  conditions  following 
Equation  (B.  lh  However,,  if  a  path  consists  of  more  than  one  subarc,  these 
constants  will  assume  different  values  because  of  the  discontinuities  at  the 
junctions  of  the  subarcs. 

The  development  assumes,  that  the  second  variation  is  evaluated  at  a  point  t3  on 
the  last  subarc.  The  assumption  of  continuity  at  this  point  requires  that 


6x(t3)  =  ni2  (t3) a  =  U(t3)b  . 


When  the  second  variation  is  evaluated  at  t  =  ,  it  is  found  that 


a»  TT12Mt3)  V(t3)  -  n22(t3)  V(t3)  b  ;>  0 


(B.  13) 


(B.  14) 


for  all  a,  b  satisfying  (B.  13),  Assuming  that  rr^2  has  an  inverse  at 
t  =  tg,  this  is  equivalent  to  stating  that  the  matrix  of 


b«  U'V  -  U>  (tt. 


b  s  0 


(B.  15) 


must  have  no  negative  eigenvalues  (although  some  can  be  zero).  Substitution 
of  (B.  8)  -  (B.  12)  into  (B.  15),  and  subsequent  rearrangement,  then  gives 
the  desired  result, 

b'  ni2^T)  nl2  nll^V  ^12  +  ^22  b  *  0  * 


<B,  16) 
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at  an  arbitrarily  selected  point  t  =  tg  on  the  last  subarc, 

* 

Notice  that  if  T  is  selected  as  tg,  the  matrix  of  <B.  16)  reduces  to  U(T)  =  Z. 
This  is  reasonable,  since  the  only  unbroken  extremal  between  the  two  end¬ 
points  is  the  original  extremal,  so  the  value  of  the  second  variation  is, 
naturally,  zero. 


! 
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0  2D  40  60  80  KX)  120  HO  !60  180  200 

TIME  -  SECONDS 


Figure  3-1.  Altitude  (h);  Velocity  (v);  Heating  Rate  (q);  Angle  of  Attack  (a) 

i  f  T  .  ,  o 

Versus  Time  for  Relative  Minimum  of  J  =  7p — —  (Q  -  q)  dT 

•  l°  A n 

where  T  =  200  sec;  t  =  0;  Q  =  195  BTU/sq  ft/sec, 
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VELOCITY,  v,  FT/  SECOND  x  10  3 


TIME -SECONDS 


Figure  3-3.  Altitude  (h);  Velocity  (v)  Heating  Rate  (q) ;  and  Angle  of  Attack  ( a ) 

1  T  T  *  *  2 

Versus  Time  for  Relative  Minimum  of  J  =  t -  I  (Q  -  q)  d  t  ; 

•  0 
T  =  550  sec;  t  =  0;  Q  =  195  BTU/sq  ft/sec. 
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HEATING  RATE  ,q,  BTU/SQ  FT/SEC 


Figure  3-4.  Altitude  (h);  Velocity  (v);  Heating  Rate  (q);  Angle  of  Attack  ( a ); 

Pilot’s  Acceleration  (a^);  and  Additional  Multiplier  (\i)  Versus 

Time  for  Non-Optimal  Re-entry  Trajectory 
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Figure  3 


-8.  Constrained  Angle  of  Attack  U  and  Unconstrained  Valued 
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Figure  4-2.  Standard  Roll  Function  <Mt) 
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Figure  4-5,  Height  for  Standard  Reference  yQ  =  -6,  4° 
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Figure  4 1 6,  Velocity  for  Standard  Reference.  yQ  =  -6,4° 
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Figure  4-7.  Re-entry  Angle  for  Standard  Reference 
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Figure  4-8.  Standard  Reference  Range 
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Figure  4-10.  Roll  Magnitude  Range  Effects 
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Figure  4-11.  Re-entry  Angle  Range  Effects 
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Figure  4-12.  Height  Range  Effects 
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Figure  4-13,  Velocity  Range  Effects 
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Figure  4-15.  Density  Change  Range  Effects 
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Figure  4-16.  4 
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Figure  4-18,  A<f>  Guidance  Sensitivities 
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Figure  4—19,  Guidance  Sensitivit.es 
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Figure  4-20.  Closed-Loop  Traj  ectory 


Figure  4-21.  Closed-Loop  Trajectory 
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Figure  4-22.  Closed-Loop  Trajectory 
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Figure  4-23.  Closed-Loop  Trajectory 
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Figure  4-25.  Closed-Loop  Trajectory 
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Figure  4-26.  Closed-Loop  Trajectory 
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Figure  4-27.  Closed-Loop  Trajectory 
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Figure  4-28.  Closed-Loop  Trajectory 
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Figure  4-29.  Closed-Loop  System  Showing  Perturbational 
Navigation  System 
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Figure  4-30.  First-Order  Navigation  System  for  a  Two-Dimensional 
Trajectory  -  Fixed  Reference 
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Figure  4-32.  Closed-Loop  Trajectory 


163 


'4  'UK** 


A  RANGE 
MILES 


TIME  AFTER  RESET  IN  SECONDS 


Figure  4-33.  Closed-Loop  Trajectory 
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Figure  4-35.  Closed-Loop  Trajectory 


Figure  4-36.  Closed-Loop  Trajectory 
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Figure  4-38.  Extremal  Reference  Trajectory 
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Figure  4-41, 


Range  and  Altitude  Sensitivity  for  Extremal 
Reference  Trajectory 


172 


ALTITUDE  ERROR  Ah  IOOO  FEET 


Figure  4-42.  Open-  and  Closed-Loop  Errors,  Altitude  Case 
Initial  Perturbation  An  =  0.  05  degree 
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Figure  4-43,  Open-  and  Closed-Loop  Errors,  Range  Case: 

Initial  Perturbation  ~  0.  05  degree 
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CONTROL  CORRECTION 


errata  SHEET  FOR  FDL-TDR-64-13,  volume  I, 

FEASIBILITY  STUDY  OF  NEW  TECHNIQUES  FOR  CONTROL  OF 
RE-ENTRY  VEHICLES 

page  8,  Equation  (2.1),  should  read 

T 

j  =  g(T,x(T))  +  J  fd(x,u)dr 

page  21  line  6  from  bottom  -  first  word  is  Equations 

page  23  line  5  from  bottom  -  first  word  is  Variables 

page  2 6  The  off-diagonal  =terms  in  matrix  at  bottom  of  page  are 
•  ♦ 

V  J  and  V'  J  respectively.. 

P  P 

ro  o 

page  27  last  equation  is  numbered  (2.71)*  Also,  the  g  in  the  last  term  should 
be  Y  . 

page  31  Equation  (2.92)  requires  brackets  as  in  equation  (2.91)* 
page  32  line  2,  (A28)  should  be  (A. 28) 

Equation  (2.96)  should  read 

(v  x  ,.)  =  V  x  (t.)  +  Cx  (t  )  -  x+  (t.)lv  t1 

v  p  n+r,  p  n+lv  1'  n+ls  V  n+lv  V  p  1 
ro  1  ro  o 

page  34  -  line  10  from  bottom  -  capital  P  should  be  lower  case  for  page 
page  35  delete  lines  9  and  10  from  bottom 

page  46  Equation  ( 3 ♦ 4l )  p  should  be  p1  (transposed,  or  row  vector). 

page  48  line  8,  a^+2  should  be 

page  58  line  8  from  bottom  -  the  a  should  be  a 

After  line  6  add  the  following:  The  minimum  principle  is  satisfied 
by  the  rales  at  the  top  of  page  49.  A  corner  exists  at  a  point  where 
Pg  passes  through  zero.  Integrations  must  be  restarted  at  such  a  point 
WlCh  <W  =  *  '  “oiJ>>  provided  that  p^  f  0. 


page  53  line  8  from  bottom  -  Ct should  be  a. 

page  59  Equation  (3.67)  -  in  the  third  equation  should  be  . 

dx3  '3 

equation  :(3*68),  last  term  -  u^  is 
page  60  line  8,  the  first,  term  of  the  vector  -is- -p^vC^  . 
page  6l  line  3j  the  3 rd  and  4th  words  should  be  "when  either" 
line  7  from  bottom  -  u^  should  be  fj 
page  63  equation  (3»7-9)  should  read  -juj  <  (f)  <  | u I 

equation  for  ij^  is  numbered  (3.86),  and  should  read 


f°L  +  CD®W  +  .  .  .. 

_ _ _ _ _  sin(u-0) 


(b^  -  ac) 


s  mu 


(3. 


line  10  should  read  ...and  -jj^  -  0  with  /j^  ^  0,.... 

Equation  (3.81)  "  the  first  column  elements  should  have  dots  above 

them  (except  for  the  zero)  to  indicate  time  derivatives,  i.e.,  v(t),  etc. 

page  64  Equation  (3.82)  -  the  first  element  of  the  matrix  is  -p/(T)x(T)  and  the 

third  element  is  -V/  x(T)p(T). 

^o 

page  65  Equation  (3*84)  The  left-hand  side  should  be  the  same  as  that  for 
equation  (2.35)  on  page  17. 

Expression  (3*85)  should  read 

(',V2)'lvV(V2>"1 

Equation  (3*86)  -  The  first  column  elements  of  the  inverted  matrix  should 
have  dots  over  them  to  indicate  time  derivations. 

Line  4  from  bottom  -  (3«85)  should  be  (3*86) 


3. 


page  66  line  2  -  "completed1’  should  be  "computed", 
page  70  -  equation  (4.7)  should  read 

Ax  =  [vf  +  v  f  E-]Ax 

x  u 

page  72  equation  (4.11)  -  the  arguments  t  in  the  second  terras  should  be  "a", 

page  73  equation  (4.14)  -  the  ao's  should  all  be  a's* 

the  dimension  of  A*  should  be  ((n-q)xn) 
page  74  -  line  6  from  bottom  -  range  on  t  should  be  a  <  t  <  T. 
page  75  equation  (4.22)  should  read  p  =  PQe~^R^ 
page  76  line  2  -  the  atmosphere  coefficient  should  be  I/235OO. 

line  4  from  bottom:.  Yq  =  -6.4  degrees. 

Pa8e  79  last  line...  the  10-g  acceleration  limit. 

page  94  line  3  -  the  last  word  is  "law" 

page  95  line  4  -  a  should  be  y 

page  110  Reference  6  -  replace  "Report"  by  TR-R- 

Reference  7  -  "Atmosphere"  should  be  "Atmospheres" 
replace  "Report"  by  TR- 
Reference  9  -  Replace  "Report"  by  TR- 
page  111  Reference  19  -  the  original  paper  number  is  6l-6 
page  112  Reference  27  -  add  ( L961) 

page  113  the  right-hand  side  of  equation  (A. l)  should  be  -p ' 

Equation  (A. 4)  should  read  0  >  pt 

page  119  the  initial  condition  for  equation  (A.36)  is  *n+1(0)  =  0 
the  lower  limit  for  the  integraL  o--c  equation  (A.37)  is  0. 
page.  122  delete  the  last  sentence  of  the  first  paragraph 


